Weierstrass elliptic function

Template:Short description Template:Redirect In mathematics, the Weierstrass elliptic functions are elliptic functions that take a particularly simple form. They are named for Karl Weierstrass. This class of functions is also referred to as ℘-functions and they are usually denoted by the symbol ℘, a uniquely fancy script p. They play an important role in the theory of elliptic functions, i.e., meromorphic functions that are doubly periodic. A ℘-function together with its derivative can be used to parameterize elliptic curves and they generate the field of elliptic functions with respect to a given period lattice.

Symbol for Weierstrass $\mathrm{\wp }$-function

Model of Weierstrass $\mathrm{\wp }$-function

Motivation

A cubic of the form $C_{g_2,g_3}^{ℂ}=\{(x,y)\in {ℂ}^2:y^2=4x^3-g_{2}x-g_3\}$, where $g_2,g_3\in ℂ$ are complex numbers with $g_2^3-27g_3^2\ne 0$, cannot be rationally parameterized. Yet one still wants to find a way to parameterize it.

For the quadric $K=\{(x,y)\in {ℝ}^2:x^2+y^2=1\}$; the unit circle, there exists a (non-rational) parameterization using the sine function and its derivative the cosine function: $$\psi :ℝ/2\pi \mathbb{Z}\to K,t\mapsto (\sin t,\cos t).$$ Because of the periodicity of the sine and cosine $ℝ/2\pi \mathbb{Z}$ is chosen to be the domain, so the function is bijective.

In a similar way one can get a parameterization of $C_{g_2,g_3}^{ℂ}$ by means of the doubly periodic $\mathrm{\wp }$-function and its derivative, namely via $(x,y)=(\mathrm{\wp }(z),{\mathrm{\wp }}^{\prime }(z))$. This parameterization has the domain $ℂ/\mathrm{\Lambda }$, which is topologically equivalent to a torus.

There is another analogy to the trigonometric functions. Consider the integral function $$a(x)={\int }_0^x\frac{dy}{\sqrt{1-y^2}}.$$ It can be simplified by substituting $y=\sin t$ and $s=\arcsin x$: $$a(x)={\int }_0^{s}dt=s=\arcsin x.$$ That means $a^{-1}(x)=\sin x$. So the sine function is an inverse function of an integral function.

Elliptic functions are the inverse functions of elliptic integrals. In particular, let: $$u(z)={\int }_z^{\mathrm{\infty }}\frac{ds}{\sqrt{4s^3-g_{2}s-g_3}}.$$ Then the extension of $u^{-1}$ to the complex plane equals the $\mathrm{\wp }$-function. This invertibility is used in complex analysis to provide a solution to certain nonlinear differential equations satisfying the Painlevé property, i.e., those equations that admit poles as their only movable singularities.

Definition

Visualization of the $\mathrm{\wp }$-function with invariants $g_2=1+i$ and $g_3=2-3i$ in which white corresponds to a pole, black to a zero.

Let ${\omega }_1,{\omega }_2\in ℂ$ be two complex numbers that are linearly independent over $ℝ$ and let $\mathrm{\Lambda }:=\mathbb{Z}{\omega }_1+\mathbb{Z}{\omega }_2:=\{m{\omega }_1+n{\omega }_2:m,n\in \mathbb{Z}\}$ be the period lattice generated by those numbers. Then the $\mathrm{\wp }$-function is defined as follows:

$\mathrm{\wp }(z,{\omega }_1,{\omega }_2):=\mathrm{\wp }(z)=\frac{1}{z^2}+\sum _{\lambda \in \mathrm{\Lambda }\setminus \{0\}}(\frac{1}{(z-\lambda {)}^2}-\frac{1}{{\lambda }^2}).$

This series converges locally uniformly absolutely in the complex torus $ℂ/\mathrm{\Lambda }$.

It is common to use $1$ and $\tau$ in the upper half-plane $ℍ:=\{z\in ℂ:\mathrm{Im}(z)>0\}$ as generators of the lattice. Dividing by ${\omega }_1$ maps the lattice $\mathbb{Z}{\omega }_1+\mathbb{Z}{\omega }_2$ isomorphically onto the lattice $\mathbb{Z}+\mathbb{Z}\tau$ with $\tau =\frac{{\omega }_2}{{\omega }_1}$. Because $-\tau$ can be substituted for $\tau$, without loss of generality we can assume $\tau \in ℍ$, and then define $\mathrm{\wp }(z,\tau ):=\mathrm{\wp }(z,1,\tau )$. With that definition, we have $\mathrm{\wp }(z,{\omega }_1,{\omega }_2)={\omega }_1^{-2}\mathrm{\wp }(z/{\omega }_1,{\omega }_2/{\omega }_1)$.

Properties

• $\mathrm{\wp }$ is a meromorphic function with a pole of order 2 at each period $\lambda$ in $\mathrm{\Lambda }$.
• $\mathrm{\wp }$ is a homogeneous function in that:

$\mathrm{\wp }(\lambda z,\lambda {\omega }_1,\lambda {\omega }_2)={\lambda }^{-2}\mathrm{\wp }(z,{\omega }_1,{\omega }_2).$

• $\mathrm{\wp }$ is an even function. That means $\mathrm{\wp }(z)=\mathrm{\wp }(-z)$ for all $z\in ℂ\setminus \mathrm{\Lambda }$, which can be seen in the following way:

$\begin{array}{cc}\mathrm{\wp }(-z)& =\frac{1}{(-z{)}^2}+\sum _{\lambda \in \mathrm{\Lambda }\setminus \{0\}}(\frac{1}{(-z-\lambda {)}^2}-\frac{1}{{\lambda }^2})\\ & =\frac{1}{z^2}+\sum _{\lambda \in \mathrm{\Lambda }\setminus \{0\}}(\frac{1}{(z+\lambda {)}^2}-\frac{1}{{\lambda }^2})\\ & =\frac{1}{z^2}+\sum _{\lambda \in \mathrm{\Lambda }\setminus \{0\}}(\frac{1}{(z-\lambda {)}^2}-\frac{1}{{\lambda }^2})=\mathrm{\wp }(z).\end{array}$

The second last equality holds because $\{-\lambda :\lambda \in \mathrm{\Lambda }\}=\mathrm{\Lambda }$. Since the sum converges absolutely this rearrangement does not change the limit.

• The derivative of $\mathrm{\wp }$ is given by: $${\mathrm{\wp }}^{\prime }(z)=-2\sum _{\lambda \in \mathrm{\Lambda }}\frac{1}{(z-\lambda {)}^3}.$$
• $\mathrm{\wp }$ and ${\mathrm{\wp }}^{\prime }$ are doubly periodic with the periods ${\omega }_1$ and ${\omega }_2$. This means: $$\begin{array}{cc}\mathrm{\wp }(z+{\omega }_1)& =\mathrm{\wp }(z)=\mathrm{\wp }(z+{\omega }_2),\text{and}\\ {\mathrm{\wp }}^{\prime }(z+{\omega }_1)& ={\mathrm{\wp }}^{\prime }(z)={\mathrm{\wp }}^{\prime }(z+{\omega }_2).\end{array}$$ It follows that $\mathrm{\wp }(z+\lambda )=\mathrm{\wp }(z)$ and ${\mathrm{\wp }}^{\prime }(z+\lambda )={\mathrm{\wp }}^{\prime }(z)$ for all $\lambda \in \mathrm{\Lambda }$.

Laurent expansion

Let $r:=min\{|\lambda |:0\ne \lambda \in \mathrm{\Lambda }\}$. Then for $0<|z|<r$ the $\mathrm{\wp }$-function has the following Laurent expansion $$\mathrm{\wp }(z)=\frac{1}{z^2}+\sum _{n=1}^{\mathrm{\infty }}(2n+1)G_{2n+2}{z}^{2n}$$ where $$G_n=\sum _{0\ne \lambda \in \mathrm{\Lambda }}{\lambda }^{-n}$$ for $n\ge 3$ are so called Eisenstein series.

Differential equation

Set $g_2=60G_4$ and $g_3=140G_6$. Then the $\mathrm{\wp }$-function satisfies the differential equation $$\mathrm{\wp }{\text{'}}^2(z)=4{\mathrm{\wp }}^3(z)-g_2\mathrm{\wp }(z)-g_3.$$ This relation can be verified by forming a linear combination of powers of $\mathrm{\wp }$ and ${\mathrm{\wp }}^{\prime }$ to eliminate the pole at $z=0$. This yields an entire elliptic function that has to be constant by Liouville's theorem.

Invariants

The real part of the invariant g₃ as a function of the square of the nome q on the unit disk.

The imaginary part of the invariant g₃ as a function of the square of the nome q on the unit disk.

The coefficients of the above differential equation $g_2$ and $g_3$ are known as the invariants. Because they depend on the lattice $\mathrm{\Lambda }$ they can be viewed as functions in ${\omega }_1$ and ${\omega }_2$.

The series expansion suggests that $g_2$ and $g_3$ are homogeneous functions of degree $-4$ and $-6$. That is $$g_2(\lambda {\omega }_1,\lambda {\omega }_2)={\lambda }^{-4}{g}_2({\omega }_1,{\omega }_2)$$ $$g_3(\lambda {\omega }_1,\lambda {\omega }_2)={\lambda }^{-6}{g}_3({\omega }_1,{\omega }_2)$$ for $\lambda \ne 0$.

If ${\omega }_1$ and ${\omega }_2$ are chosen in such a way that $\mathrm{Im}\left(\frac{{\omega }_2}{{\omega }_1}\right)>0$, $g_2$ and $g_3$ can be interpreted as functions on the upper half-plane $ℍ:=\{z\in ℂ:\mathrm{Im}(z)>0\}$.

Let $\tau =\frac{{\omega }_2}{{\omega }_1}$. One has: $$g_2(1,\tau )={\omega }_1^4{g}_2({\omega }_1,{\omega }_2),$$ $$g_3(1,\tau )={\omega }_1^6{g}_3({\omega }_1,{\omega }_2).$$ That means g₂ and g₃ are only scaled by doing this. Set $$g_2(\tau ):=g_2(1,\tau )$$ and $$g_3(\tau ):=g_3(1,\tau ).$$ As functions of $\tau \in ℍ$, $g_2$ and $g_3$ are so called modular forms.

The Fourier series for $g_2$ and $g_3$ are given as follows: $$g_2(\tau )=\frac{4}{3}{\pi }^4[1+240\sum _{k=1}^{\mathrm{\infty }}{\sigma }_3(k)q^{2k}]$$ $$g_3(\tau )=\frac{8}{27}{\pi }^6[1-504\sum _{k=1}^{\mathrm{\infty }}{\sigma }_5(k)q^{2k}]$$ where $${\sigma }_m(k):=\sum _{d\mid k}{d}^m$$ is the divisor function and $q=e^{\pi i\tau }$ is the nome.

Modular discriminant

The real part of the discriminant as a function of the square of the nome q on the unit disk.

The modular discriminant $\mathrm{\Delta }$ is defined as the discriminant of the characteristic polynomial of the differential equation $$\mathrm{\wp }{\text{'}}^2(z)=4{\mathrm{\wp }}^3(z)-g_2\mathrm{\wp }(z)-g_3$$ as follows: $$\mathrm{\Delta }=g_2^3-27g_3^2.$$ The discriminant is a modular form of weight $12$. That is, under the action of the modular group, it transforms as $$\mathrm{\Delta }\left(\frac{a\tau +b}{c\tau +d}\right)={(c\tau +d)}^{12}\mathrm{\Delta }(\tau )$$ where $a,b,d,c\in \mathbb{Z}$ with $ad-bc=1$.

Note that $\mathrm{\Delta }=(2\pi {)}^{12}{\eta }^{24}$ where $\eta$ is the Dedekind eta function.

For the Fourier coefficients of $\mathrm{\Delta }$, see Ramanujan tau function.

The constants e₁, e₂ and e₃

$e_1$, $e_2$ and $e_3$ are usually used to denote the values of the $\mathrm{\wp }$-function at the half-periods. $$e_1\equiv \mathrm{\wp }\left(\frac{{\omega }_1}{2}\right)$$ $$e_2\equiv \mathrm{\wp }\left(\frac{{\omega }_2}{2}\right)$$ $$e_3\equiv \mathrm{\wp }\left(\frac{{\omega }_1+{\omega }_2}{2}\right)$$ They are pairwise distinct and only depend on the lattice $\mathrm{\Lambda }$ and not on its generators.

$e_1$, $e_2$ and $e_3$ are the roots of the cubic polynomial $4\mathrm{\wp }(z{)}^3-g_2\mathrm{\wp }(z)-g_3$ and are related by the equation: $$e_1+e_2+e_3=0.$$ Because those roots are distinct the discriminant $\mathrm{\Delta }$ does not vanish on the upper half plane. Now we can rewrite the differential equation: $$\mathrm{\wp }{\text{'}}^2(z)=4(\mathrm{\wp }(z)-e_1)(\mathrm{\wp }(z)-e_2)(\mathrm{\wp }(z)-e_3).$$ That means the half-periods are zeros of ${\mathrm{\wp }}^{\prime }$.

The invariants $g_2$ and $g_3$ can be expressed in terms of these constants in the following way: $$g_2=-4(e_1{e}_2+e_1{e}_3+e_2{e}_3)$$ $$g_3=4e_1{e}_2{e}_3$$ $e_1$, $e_2$ and $e_3$ are related to the modular lambda function: $$\lambda (\tau )=\frac{e_3-e_2}{e_1-e_2},\tau =\frac{{\omega }_2}{{\omega }_1}.$$

Relation to Jacobi's elliptic functions

For numerical work, it is often convenient to calculate the Weierstrass elliptic function in terms of Jacobi's elliptic functions.

The basic relations are: $$\mathrm{\wp }(z)=e_3+\frac{e_1-e_3}{{\mathrm{sn}}^{2}w}=e_2+(e_1-e_3)\frac{{\mathrm{dn}}^{2}w}{{\mathrm{sn}}^{2}w}=e_1+(e_1-e_3)\frac{{\mathrm{cn}}^{2}w}{{\mathrm{sn}}^{2}w}$$ where $e_1,e_2$ and $e_3$ are the three roots described above and where the modulus k of the Jacobi functions equals $$k=\sqrt{\frac{e_2-e_3}{e_1-e_3}}$$ and their argument w equals $$w=z\sqrt{e_1-e_3}.$$

Relation to Jacobi's theta functions

The function $\mathrm{\wp }(z,\tau )=\mathrm{\wp }(z,1,{\omega }_2/{\omega }_1)$ can be represented by Jacobi's theta functions: $$\mathrm{\wp }(z,\tau )={(\pi {\theta }_2(0,q){\theta }_3(0,q)\frac{{\theta }_4(\pi z,q)}{{\theta }_1(\pi z,q)})}^2-\frac{{\pi }^2}{3}({\theta }_2^4(0,q)+{\theta }_3^4(0,q))$$ where $q=e^{\pi i\tau }$ is the nome and $\tau$ is the period ratio $(\tau \in ℍ)$. This also provides a very rapid algorithm for computing $\mathrm{\wp }(z,\tau )$.

Relation to elliptic curves

Template:See also Consider the embedding of the cubic curve in the complex projective plane

${\overline{C}}_{g_2,g_3}^{ℂ}=\{(x,y)\in {ℂ}^2:y^2=4x^3-g_{2}x-g_3\}\cup \{O\}\subset {ℂ}^2\cup {\mathbb{P}}_1(ℂ)={\mathbb{P}}_2(ℂ).$

where $O$ is a point lying on the line at infinity ${\mathbb{P}}_1(ℂ)$. For this cubic there exists no rational parameterization, if $\mathrm{\Delta }\ne 0$. In this case it is also called an elliptic curve. Nevertheless there is a parameterization in homogeneous coordinates that uses the $\mathrm{\wp }$-function and its derivative ${\mathrm{\wp }}^{\prime }$:

$\phi (\mathrm{\wp },{\mathrm{\wp }}^{\prime }):ℂ/\mathrm{\Lambda }\to {\overline{C}}_{g_2,g_3}^{ℂ},z\mapsto \{\begin{array}{cc}[\mathrm{\wp }(z):{\mathrm{\wp }}^{\prime }(z):1]& z\notin \mathrm{\Lambda }\\ [0:1:0]& z\in \mathrm{\Lambda }\end{array}$

Now the map $\phi$ is bijective and parameterizes the elliptic curve ${\overline{C}}_{g_2,g_3}^{ℂ}$.

$ℂ/\mathrm{\Lambda }$ is an abelian group and a topological space, equipped with the quotient topology.

It can be shown that every Weierstrass cubic is given in such a way. That is to say that for every pair $g_2,g_3\in ℂ$ with $\mathrm{\Delta }=g_2^3-27g_3^2\ne 0$ there exists a lattice $\mathbb{Z}{\omega }_1+\mathbb{Z}{\omega }_2$, such that

$g_2=g_2({\omega }_1,{\omega }_2)$ and $g_3=g_3({\omega }_1,{\omega }_2)$.

Addition theorem

The addition theorem states that if $z,w,$ and $z+w$ do not belong to $\mathrm{\Lambda }$, then $$det\left[\begin{array}{ccc}1& \mathrm{\wp }(z)& {\mathrm{\wp }}^{\prime }(z)\\ 1& \mathrm{\wp }(w)& {\mathrm{\wp }}^{\prime }(w)\\ 1& \mathrm{\wp }(z+w)& -{\mathrm{\wp }}^{\prime }(z+w)\end{array}\right]=0.$$ This states that the points $P=(\mathrm{\wp }(z),{\mathrm{\wp }}^{\prime }(z)),$ $Q=(\mathrm{\wp }(w),{\mathrm{\wp }}^{\prime }(w)),$ and $R=(\mathrm{\wp }(z+w),-{\mathrm{\wp }}^{\prime }(z+w))$ are collinear, the geometric form of the group law of an elliptic curve.

This can be proven by considering constants $A,B$ such that $${\mathrm{\wp }}^{\prime }(z)=A\mathrm{\wp }(z)+B,{\mathrm{\wp }}^{\prime }(w)=A\mathrm{\wp }(w)+B.$$ Then the elliptic function $${\mathrm{\wp }}^{\prime }(\zeta )-A\mathrm{\wp }(\zeta )-B$$ has a pole of order three at zero, and therefore three zeros whose sum belongs to $\mathrm{\Lambda }$. Two of the zeros are $z$ and $w$, and thus the third is congruent to $-z-w$.

Alternative form

The addition theorem can be put into the alternative form, for $z,w,z-w,z+w\in ̸\mathrm{\Lambda }$: $$\mathrm{\wp }(z+w)=\frac{1}{4}{\left[\frac{{\mathrm{\wp }}^{\prime }(z)-{\mathrm{\wp }}^{\prime }(w)}{\mathrm{\wp }(z)-\mathrm{\wp }(w)}\right]}^2-\mathrm{\wp }(z)-\mathrm{\wp }(w).$$

As well as the duplication formula: $$\mathrm{\wp }(2z)=\frac{1}{4}{\left[\frac{{\mathrm{\wp }}^{″}(z)}{{\mathrm{\wp }}^{\prime }(z)}\right]}^2-2\mathrm{\wp }(z).$$

Proofs

This can be proven from the addition theorem shown above. The points $P=(\mathrm{\wp }(u),{\mathrm{\wp }}^{\prime }(u)),Q=(\mathrm{\wp }(v),{\mathrm{\wp }}^{\prime }(v)),$ and $R=(\mathrm{\wp }(u+v),-{\mathrm{\wp }}^{\prime }(u+v))$ are collinear and lie on the curve $y^2=4x^3-g_{2}x-g_3$. The slope of that line is $$m=\frac{y_P-y_Q}{x_P-x_Q}=\frac{{\mathrm{\wp }}^{\prime }(u)-{\mathrm{\wp }}^{\prime }(v)}{\mathrm{\wp }(u)-\mathrm{\wp }(v)}.$$ So $x=x_P=\mathrm{\wp }(u)$, $x=x_Q=\mathrm{\wp }(v)$, and $x=x_R=\mathrm{\wp }(u+v)$ all satisfy a cubic $$(mx+q{)}^2=4x^3-g_{2}x-g_3,$$ where $q$ is a constant. This becomes $$4x^3-m^2{x}^2-(2mq+g_2)x-g_3-q^2=0.$$ Thus $x_P+x_Q+x_R=\frac{m^2}{4}$ which provides the wanted formula $\mathrm{\wp }(u+v)+\mathrm{\wp }(u)+\mathrm{\wp }(v)=\frac{1}{4}{\left[\frac{{\mathrm{\wp }}^{\prime }(u)-{\mathrm{\wp }}^{\prime }(v)}{\mathrm{\wp }(u)-\mathrm{\wp }(v)}\right]}^2.$

A direct proof is as follows. Any elliptic function $f$ can be expressed as: $$f(u)=c\prod _{i=1}^n\frac{\sigma (u-a_i)}{\sigma (u-b_i)}c\in ℂ$$ where $\sigma$ is the Weierstrass sigma function and $a_i,b_i$ are the respective zeros and poles in the period parallelogram. Considering the function $\mathrm{\wp }(u)-\mathrm{\wp }(v)$ as a function of $u$, we have $$\mathrm{\wp }(u)-\mathrm{\wp }(v)=c\frac{\sigma (u+v)\sigma (u-v)}{\sigma (u{)}^2}.$$ Multiplying both sides by $u^2$ and letting $u\to 0$, we have $1=-c\sigma (v{)}^2$, so $c=-\frac{1}{\sigma (v{)}^2}⟹\mathrm{\wp }(u)-\mathrm{\wp }(v)=-\frac{\sigma (u+v)\sigma (u-v)}{\sigma (u{)}^2\sigma (v{)}^2}.$

By definition the Weierstrass zeta function: $\frac{d}{dz}\ln \sigma (z)=\zeta (z)$ therefore we logarithmically differentiate both sides with respect to $u$ obtaining: $$\frac{{\mathrm{\wp }}^{\prime }(u)}{\mathrm{\wp }(u)-\mathrm{\wp }(v)}=\zeta (u+v)-2\zeta (u)-\zeta (u-v)$$ Once again by definition ${\zeta }^{\prime }(z)=-\mathrm{\wp }(z)$ thus by differentiating once more on both sides and rearranging the terms we obtain $$-\mathrm{\wp }(u+v)=-\mathrm{\wp }(u)+\frac{1}{2}\frac{{\mathrm{\wp }}^{″}(v)[\mathrm{\wp }(u)-\mathrm{\wp }(v)]-{\mathrm{\wp }}^{\prime }(u)[{\mathrm{\wp }}^{\prime }(u)-{\mathrm{\wp }}^{\prime }(v)]}{[\mathrm{\wp }(u)-\mathrm{\wp }(v){]}^2}$$ Knowing that ${\mathrm{\wp }}^{″}$ has the following differential equation $2{\mathrm{\wp }}^{″}=12{\mathrm{\wp }}^2-g_2$ and rearranging the terms one gets the wanted formula $$\mathrm{\wp }(u+v)=\frac{1}{4}{\left[\frac{{\mathrm{\wp }}^{\prime }(u)-{\mathrm{\wp }}^{\prime }(v)}{\mathrm{\wp }(u)-\mathrm{\wp }(v)}\right]}^2-\mathrm{\wp }(u)-\mathrm{\wp }(v).$$

Typography

The Weierstrass's elliptic function is usually written with a rather special, lower case script letter ℘, which was Weierstrass's own notation introduced in his lectures of 1862–1863.Template:Refn It should not be confused with the normal mathematical script letters P: 𝒫 and 𝓅.

In computing, the letter ℘ is available as \wp in TeX. In Unicode the code point is Template:Unichar, with the more correct alias Template:Smallcaps.Template:Refn In HTML, it can be escaped as ℘ or ℘. Template:Charmap

See also

• Weierstrass functions
• Jacobi elliptic functions
• Lemniscate elliptic functions

References

• Template:AS ref
• N. I. Akhiezer, Elements of the Theory of Elliptic Functions, (1970) Moscow, translated into English as AMS Translations of Mathematical Monographs Volume 79 (1990) AMS, Rhode Island Template:Isbn
• Tom M. Apostol, Modular Functions and Dirichlet Series in Number Theory, Second Edition (1990), Springer, New York Template:Isbn (See chapter 1.)
• K. Chandrasekharan, Elliptic functions (1980), Springer-Verlag Template:Isbn
• Konrad Knopp, Funktionentheorie II (1947), Dover Publications; Republished in English translation as Theory of Functions (1996), Dover Publications Template:Isbn
• Serge Lang, Elliptic Functions (1973), Addison-Wesley, Template:Isbn
• E. T. Whittaker and G. N. Watson, A Course of Modern Analysis, Cambridge University Press, 1952, chapters 20 and 21

External links

• Template:Springer
• Weierstrass's elliptic functions on Mathworld.
• Chapter 23, Weierstrass Elliptic and Modular Functions in DLMF (Digital Library of Mathematical Functions) by W. P. Reinhardt and P. L. Walker.
• Weierstrass P function and its derivative implemented in C by David Dumas