Pi, tau auxillary functions. still bugs in the poles.

Former-commit-id: 99d550f5ffa03d001f89e41a07e70652105bdc36

notes/Scattering and Shit.lyx

@@ -240,6 +240,117 @@ Orthonormality:
 Pi and tau
 \end_layout
 
+\begin_layout Subsection
+Xu
+\begin_inset CommandInset label
+LatexCommand label
+name "sub:Xu pitau"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+As in (37)
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula 
+\begin{eqnarray*}
+\pi_{mn}\left(\cos\theta\right) & = & \frac{m}{\sin\theta}P_{n}^{m}\left(\cos\theta\right)\\
+\tau_{mn}\left(\cos\theta\right) & = & \frac{\ud}{\ud\theta}P_{n}^{m}\left(\cos\theta\right)=-\left(\sin\theta\right)\frac{\ud P_{n}^{m}\left(\cos\theta\right)}{\ud\left(\cos\theta\right)}
+\end{eqnarray*}
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+The expressions 
+\begin_inset Formula $\left(\sin\theta\right)^{-1}$
+\end_inset
+
+ and 
+\begin_inset Formula $\frac{\ud P_{n}^{m}\left(\cos\theta\right)}{\ud\left(\cos\theta\right)}$
+\end_inset
+
+ are singular for 
+\begin_inset Formula $\cos\theta=\pm1$
+\end_inset
+
+, the limits 
+\begin_inset Formula $\tau_{mn}\left(\pm1\right),\pi_{mn}\left(\pm1\right)$
+\end_inset
+
+ however exist.
+ Labeling 
+\begin_inset Formula $x\equiv\cos\theta$
+\end_inset
+
+, 
+\begin_inset Formula $\sqrt{\left(1+x\right)\left(1-x\right)}=\sqrt{1-x^{2}}\equiv\sin\theta$
+\end_inset
+
+ and using the asymptotic expression (DLMF 14.8.2) we obtain that the limits
+ are nonzero only for 
+\begin_inset Formula $m=\pm1$
+\end_inset
+
+ and
+\begin_inset Formula 
+\begin{eqnarray*}
+\pi_{1\nu}(1-) & = & -\frac{\nu\left(\nu+1\right)\left(\nu+2\right)}{2}\\
+\tau_{1\nu}(1-) & = & \frac{\nu\left(\nu+1\right)\left(\nu+2\right)}{2}
+\end{eqnarray*}
+
+\end_inset
+
+and using the parity property 
+\begin_inset Formula $P_{n}^{m}\left(-x\right)=\left(-1\right)^{m+n}P_{n}^{m}\left(x\right)$
+\end_inset
+
+
+\begin_inset Formula 
+\begin{eqnarray*}
+\pi_{1\nu}(-1+) & = & \left(-1\right)^{\nu}\frac{\nu\left(\nu+1\right)\left(\nu+2\right)}{2}\\
+\tau_{1\nu}(-1+) & = & \left(-1\right)^{\nu}\frac{\nu\left(\nu+1\right)\left(\nu+2\right)}{2}
+\end{eqnarray*}
+
+\end_inset
+
+For 
+\begin_inset Formula $m=1$
+\end_inset
+
+, we simply use the relation 
+\begin_inset Formula $P_{n}^{-m}=\left(CS\right)^{m}P_{n}^{m}\frac{\left(n-m\right)!}{\left(n+m\right)!}$
+\end_inset
+
+ to get
+\begin_inset Formula 
+\begin{eqnarray*}
+\pi_{-1\nu}(1-) & = & -CS\frac{\nu+2}{2}\\
+\tau_{-1\nu}(1-) & = & CS\frac{\nu+2}{2}\\
+\pi_{-1\nu}(-1+) & = & CS\left(-1\right)^{\nu}\frac{\nu+2}{2}\\
+\tau_{-1\nu}(-1+) & = & CS\left(-1\right)^{\nu}\frac{\nu+2}{2}
+\end{eqnarray*}
+
+\end_inset
+
+where 
+\begin_inset Formula $CS$
+\end_inset
+
+ is 
+\begin_inset Formula $-1$
+\end_inset
+
+ if the Condon-Shortley phase is employed on the level of Legendre polynomials,
+ 1 otherwise.
+\end_layout
+
 \begin_layout Subsection
 Taylor
 \end_layout
@@ -256,6 +367,51 @@ Taylor
 
 \end_layout
 
+\begin_layout Standard
+The limiting expressions are obtained simply by multiplying the expressions
+ from sec.
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sub:Xu pitau"
+
+\end_inset
+
+ by the normalisation factor,
+\begin_inset Formula 
+\begin{eqnarray*}
+\tilde{\pi}_{1\nu}(1-) & = & -\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}\\
+\tilde{\tau}_{1\nu}(1-) & = & \sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}\\
+\tilde{\pi}_{1\nu}(-1+) & = & \left(-1\right)^{\nu}\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}\\
+\tilde{\tau}_{1\nu}(-1+) & = & \left(-1\right)^{\nu}\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}
+\end{eqnarray*}
+
+\end_inset
+
+
+\begin_inset Formula 
+\begin{eqnarray*}
+\tilde{\pi}_{-1\nu}(1-) & = & -CS\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}\\
+\tilde{\tau}_{-1\nu}(1-) & = & CS\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}\\
+\tilde{\pi}_{-1\nu}(-1+) & = & CS\left(-1\right)^{\nu}\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}\\
+\tilde{\tau}_{-1\nu}(-1+) & = & CS\left(-1\right)^{\nu}\sqrt{\frac{2\nu+1}{4\pi}}\frac{\sqrt{\nu\left(\nu+1\right)}\left(\nu+2\right)}{2}
+\end{eqnarray*}
+
+\end_inset
+
+i.e.
+ the expressions for 
+\begin_inset Formula $m=-1$
+\end_inset
+
+ are the same as for 
+\begin_inset Formula $m=1$
+\end_inset
+
+ except for the sign if Condon-Shortley phase is used on the Legendre polynomial
+ level.
+\end_layout
+
 \begin_layout Section
 Vector spherical harmonics (?)
 \end_layout
@@ -550,6 +706,27 @@ Definition [T](2.40);
 \end_inset
 
 
+\end_layout
+
+\begin_layout Subsection
+Xu
+\end_layout
+
+\begin_layout Standard
+are the electric and magnetic waves, respectively:
+\end_layout
+
+\begin_layout Standard
+\begin_inset Formula 
+\begin{eqnarray*}
+\vect N_{mn}^{(j)} & = & \frac{n(n+1)}{kr}P_{n}^{m}\left(\cos\theta\right)e^{im\phi}z_{n}^{j}\left(kr\right)\hat{\vect r}\\
+ &  & +\left[\tau_{mn}\left(\cos\theta\right)\hat{\vect{\theta}}+i\pi_{mn}\left(\cos\theta\right)\hat{\vect{\phi}}\right]e^{im\phi}\frac{1}{kr}\frac{\ud\left(kr\, z_{n}^{j}\left(kr\right)\right)}{\ud(kr)}\\
+\vect M_{mn}^{(j)} & = & \left[i\pi_{mn}\left(\cos\theta\right)\hat{\vect{\theta}}-\tau_{mn}\left(\cos\theta\right)\hat{\vect{\phi}}\right]e^{im\phi}z_{n}^{j}\left(kr\right)
+\end{eqnarray*}
+
+\end_inset
+
+
 \end_layout
 
 \begin_layout Subsection

qpms/indexing.h

@@ -2,9 +2,10 @@
 #define QPMS_INDEXING_H
 
 #include "qpms_types.h"
+#include <math.h>
 
 static inline qpms_y_t qpms_mn2y(qpms_m_t m, qpms_l_t n) {
-        return (qpms_y_t) n * (n + 1) + m - 1;
+        return n * (n + 1) + m - 1;
 }
 
 static inline qpms_lm_t qpms_y2n(qpms_y_t y) {

qpms/qpms_types.h

@@ -7,29 +7,25 @@
 
 // integer index types
 typedef int qpms_lm_t;
-typedef unsigned int qpms_l_t;
+typedef int qpms_l_t; // can't be unsigned because of the behaviour under - operator
 typedef qpms_lm_t qpms_m_t;
 typedef size_t qpms_y_t;
 
 typedef enum {
-	QPMS_SUCCESS = 0;
+	QPMS_SUCCESS = 0,
-	QPMS_ERROR = 1;
+	QPMS_ERROR = 1
 } qpms_errno_t;
 
 // Normalisations
 typedef enum {
 	// As in TODO
-	QPMS_NORMALIZATION_XU = 3, // NI!
+	QPMS_NORMALISATION_XU = 3, // NI!
 	// As in http://www.eit.lth.se/fileadmin/eit/courses/eit080f/Literature/book.pdf, power-normalised
-	QPMS_NORMALIZATION_KRISTENSSON = 2, // NI!
+	QPMS_NORMALISATION_KRISTENSSON = 2, // NI!
-	QPMS_NORMALIZATION_POWER = QPMS_NORMALIZATION_KRISTENSSON, // NI!
+	QPMS_NORMALISATION_POWER = QPMS_NORMALISATION_KRISTENSSON, // NI!
-	QPMS_NORMALIZATION_TAYLOR = 1,
+	QPMS_NORMALISATION_TAYLOR = 1,
-	QPMS_NORMALIZATION_UNDEF = 0
+	QPMS_NORMALISATION_UNDEF = 0
-} qpms_normalization_t;
+} qpms_normalisation_t;
-
-typedef enum {
-	QPMS_SUCCESS = 0;
-} qpms_errno_t;
 
 typedef enum {
         QPMS_BESSEL_REGULAR = 1, // regular function j

qpms/test_vswf.c

@@ -0,0 +1,25 @@
+#include "vswf.h"
+#include "indexing.h"
+#include <stdio.h>
+#include <gsl/gsl_math.h>
+const double dtheta = 0.01 * M_PI;
+const qpms_l_t lMax = 3;
+
+int main() {
+	qpms_y_t nelem = qpms_lMax2nelem(lMax);
+	for (double theta = 0.; theta <= M_PI; theta += dtheta) {
+		printf("%.5e ", theta);
+		for(qpms_normalisation_t norm = 1; norm <= 3; ++norm) {//fujka :D
+			qpms_pitau_t pt = qpms_pitau_get(theta, lMax, norm);
+			for (qpms_y_t y = 0; y < nelem; ++y)
+				printf("%.5e %.5e %.5e ", pt.leg[y], pt.pi[y], pt.tau[y]);
+			qpms_pitau_free(pt);
+		}
+		printf("\n");
+	}
+}
+
+
+
+
+

qpms/vswf.c

@@ -1,28 +1,37 @@
 #include "vswf.h"
+#include "indexing.h"
 #include <math.h>
+#include <gsl/gsl_math.h>
+#include <gsl/gsl_sf_legendre.h>
+#include <stdlib.h>
+#include <string.h>
+
+#ifndef CSPHASE
+#define CSPHASE (1.) // FIXME this should be later determined by qpms_normalization_t
+#endif
 
 // Legendre functions also for negative m, see DLMF 14.9.3
-int qpms_errno_t qpms_legendre_deriv_y_fill(double *target, double *target_deriv, double x, qpms_l_t lMax,
+qpms_errno_t qpms_legendre_deriv_y_fill(double *target, double *target_deriv, double x, qpms_l_t lMax,
 		gsl_sf_legendre_t lnorm, double csphase)
 {
-	size_t n = gsl_sf_legenre_array_n(lMax);
+	size_t n = gsl_sf_legendre_array_n(lMax);
 	double *legendre_tmp = malloc(n * sizeof(double));
 	double *legendre_deriv_tmp = malloc(n * sizeof(double));
 	int gsl_errno = gsl_sf_legendre_deriv_array_e(
-			lnorm, (size_t)lMax, x, csphase, legendre_tmp,legendre_tmp_deriv);
+			lnorm, (size_t)lMax, x, csphase, legendre_tmp,legendre_deriv_tmp);
-	for (qpms_l_t l = 0; l <= lMax; ++l)
+	for (qpms_l_t l = 1; l <= lMax; ++l)
 			for (qpms_m_t m = 0; m <= l; ++m) {
 				qpms_y_t y = qpms_mn2y(m,l);
-				size_t i = gsl_sf_legenre_array_index(l,m);
+				size_t i = gsl_sf_legendre_array_index(l,m);
 				target[y] = legendre_tmp[i];
 				target_deriv[y] = legendre_deriv_tmp[i];
 			}
 	switch(lnorm) {
-		case GSL_SF_LEGEDRE_NONE:
+		case GSL_SF_LEGENDRE_NONE:
-			for (qpms_l_t l = 0; l <= lMax; ++l)
+			for (qpms_l_t l = 1; l <= lMax; ++l)
 				for (qpms_m_t m = 1; m <= l; ++m) {
 				qpms_y_t y = qpms_mn2y(-m,l);
-				size_t i = gsl_sf_legenre_array_index(l,m);
+				size_t i = gsl_sf_legendre_array_index(l,m);
 				// viz DLMF 14.9.3, čert ví, jak je to s cs fasí.
 				double factor = exp(lgamma(l-m+1)-lgamma(n+m+1))*((m%2)?-1:1);
 				target[y] = factor * legendre_tmp[i];
@@ -32,10 +41,10 @@ int qpms_errno_t qpms_legendre_deriv_y_fill(double *target, double *target_deriv
 		case GSL_SF_LEGENDRE_SCHMIDT:
 		case GSL_SF_LEGENDRE_SPHARM:
 		case GSL_SF_LEGENDRE_FULL:
-			for (qpms_l_t l = 0; l <= lMax; ++l)
+			for (qpms_l_t l = 1; l <= lMax; ++l)
 				for (qpms_m_t m = 1; m <= l; ++m) {
 				qpms_y_t y = qpms_mn2y(-m,l);
-				size_t i = gsl_sf_legenre_array_index(l,m);
+				size_t i = gsl_sf_legendre_array_index(l,m);
 				// viz DLMF 14.9.3, čert ví, jak je to s cs fasí.
 				double factor = ((m%2)?-1:1); // this is the difference from the unnormalised case
 				target[y] = factor * legendre_tmp[i];
@@ -51,13 +60,101 @@ int qpms_errno_t qpms_legendre_deriv_y_fill(double *target, double *target_deriv
 	return QPMS_SUCCESS;
 }
 
-int qpms_legendre_deriv_y_get(double **target, double **dtarget, double x, qpms_l_t lMax, gsl_sf_legendre_t lnorm,
+qpms_errno_t qpms_legendre_deriv_y_get(double **target, double **dtarget, double x, qpms_l_t lMax, gsl_sf_legendre_t lnorm,
 		double csphase)
 {
 
 	*target = malloc(sizeof(double)*qpms_lMax2nelem(lMax));
 	*dtarget = malloc(sizeof(double)*qpms_lMax2nelem(lMax));
-	return qpms_legendre_deriv_y_fill(ar, x, lMax, lnorm, csphase);
+	return qpms_legendre_deriv_y_fill(*target, *dtarget, x, lMax, lnorm, csphase);
 }
 
 
+qpms_pitau_t qpms_pitau_get(double theta, qpms_l_t lMax, qpms_normalisation_t norm) 
+{
+	qpms_pitau_t res;
+	qpms_y_t nelem = qpms_lMax2nelem(lMax);
+	res.pi = malloc(nelem * sizeof(double));
+	res.tau = malloc(nelem * sizeof(double));
+	double ct = cos(theta), st = sin(theta);
+	if (1 == fabs(ct)) { // singular case, use DLMF 14.8.2
+		memset(res.pi, 0, nelem*sizeof(double));
+		memset(res.tau, 0, nelem*sizeof(double));
+		res.leg = calloc(nelem, sizeof(double));
+		switch(norm) {
+			case QPMS_NORMALISATION_XU:
+				for (qpms_l_t l = 1; l <= lMax; ++l) {
+					res.leg[qpms_mn2y(0, l)] = (l%2)?ct:1.;
+					double p = l*(l+1)*(l+2)/2;
+					double n = (l+2)/2.;
+					int lpar = (l%2)?-1:1;
+					res.pi [qpms_mn2y(+1, l)] = ((ct>0) ? -1 : lpar) * p;
+					res.pi [qpms_mn2y(-1, l)] = ((ct>0) ? -1 : lpar) * n * CSPHASE;
+					res.tau[qpms_mn2y(+1, l)] = ((ct>0) ? +1 : lpar) * p;
+					res.tau[qpms_mn2y(-1, l)] = ((ct>0) ? +1 : lpar) * n * CSPHASE;
+				}
+				break;
+			case QPMS_NORMALISATION_TAYLOR:
+				for (qpms_l_t l = 1; l <= lMax; ++l) {
+					res.leg[qpms_mn2y(0, l)] = ((l%2)?ct:1.)*sqrt((2*l+1)*0.25*M_1_PI);
+					int lpar = (l%2)?-1:1;
+					double fl = 0.25 * (l+2) * sqrt((2*l+1)*l*(l+1)*M_1_PI);
+					res.pi [qpms_mn2y(+1, l)] = ((ct>0) ? -1 : lpar) * fl;
+					res.pi [qpms_mn2y(-1, l)] = ((ct>0) ? -1 : lpar) * fl * CSPHASE;
+					res.tau[qpms_mn2y(+1, l)] = ((ct>0) ? +1 : lpar) * fl;
+					res.tau[qpms_mn2y(-1, l)] = ((ct>0) ? +1 : lpar) * fl * CSPHASE;
+				}
+				break;
+			case QPMS_NORMALISATION_POWER:
+				for (qpms_l_t l = 1; l <= lMax; ++l) {
+					res.leg[qpms_mn2y(0, l)] = ((l%2)?ct:1.)*sqrt((2*l+1)/(4*M_PI *l*(l+1)));
+					int lpar = (l%2)?-1:1;
+					double fl = 0.25 * (l+2) * sqrt((2*l+1)*M_1_PI);
+					res.pi [qpms_mn2y(+1, l)] = ((ct>0) ? -1 : lpar) * fl;
+					res.pi [qpms_mn2y(-1, l)] = ((ct>0) ? -1 : lpar) * fl * CSPHASE;
+					res.tau[qpms_mn2y(+1, l)] = ((ct>0) ? +1 : lpar) * fl;
+					res.tau[qpms_mn2y(-1, l)] = ((ct>0) ? +1 : lpar) * fl * CSPHASE;
+	
+				}
+				break;
+			default:
+				abort();
+		}
+	}
+	else { // cos(theta) in (-1,1), use normal calculation
+		double *legder = malloc(sizeof(double)*qpms_lMax2nelem(lMax));
+		res.leg = malloc(sizeof(double)*qpms_lMax2nelem(lMax));
+		if (qpms_legendre_deriv_y_fill(res.leg, legder, ct, lMax, 
+					norm == QPMS_NORMALISATION_XU ? GSL_SF_LEGENDRE_NONE
+						: GSL_SF_LEGENDRE_SPHARM, CSPHASE)) 
+			abort();
+		if (norm == QPMS_NORMALISATION_POWER) 
+			/* for Xu (=non-normalized) and Taylor (=sph. harm. normalized) 
+			 * the correct normalisation is already obtained from gsl_sf_legendre_deriv_array_e().
+			 * However, Kristensson ("power") normalisation differs from Taylor
+			 * by 1/sqrt(l*(l+1)) factor.
+			 */
+			for (qpms_l_t l = 1; l <= lMax; ++l) {
+				double prefac = 1./sqrt(l*(l+1));
+				for (qpms_m_t m = -l; m <= l; ++m) {
+					res.leg[qpms_mn2y(m,l)] *= prefac;
+					legder[qpms_mn2y(m,l)] *= prefac;
+				}
+			}
+		for (qpms_l_t l = 1; l <= lMax; ++l) {
+			for (qpms_m_t m = -l; m <= l; ++m) {
+						res.pi [qpms_mn2y(m,l)] = m / st * res.leg[qpms_mn2y(m,l)];
+						res.tau[qpms_mn2y(m,l)] = - st * legder[qpms_mn2y(m,l)];
+			}
+		}
+		free(legder);
+	}
+	res.lMax = lMax;
+	return res;
+}
+
+void qpms_pitau_free(qpms_pitau_t x) {
+	free(x.leg);
+	free(x.pi);
+	free(x.tau);
+}

qpms/vswf.h

@@ -5,10 +5,10 @@
 
 // Electric wave N; NI
 csphvec_t qpms_vswf_single_el(int m, int n, sph_t kdlj,
-                qpms_bessel_t btyp, qpms_normalization_t norm);
+                qpms_bessel_t btyp, qpms_normalisation_t norm);
 // Magnetic wave M; NI
 csphvec_t qpms_vswf_single_mg(int m, int n, sph_t kdlj,
-                qpms_bessel_t btyp, qpms_normalization_t norm);
+                qpms_bessel_t btyp, qpms_normalisation_t norm);
 
 // Set of electric and magnetic VSWF in spherical coordinate basis
 typedef struct {
@@ -26,29 +26,31 @@ typedef struct {
  * ( gsl/specfunc/legendre_source.c and 7.24.2 of gsl docs
  */
 
-int qpms_legendre_deriv_y_get(double **result, double **result_deriv, double x, qpms_l_t lMax, 
+qpms_errno_t qpms_legendre_deriv_y_get(double **result, double **result_deriv, double x, qpms_l_t lMax, 
-		gsle_sf_legendre_t lnorm, double csphase);
+		gsl_sf_legendre_t lnorm, double csphase); // free() result and result_deriv yourself!
-int qpms_errno_t qpms_legendre_deriv_y_fill(double *where, double *where_deriv, double x, 
+qpms_errno_t qpms_legendre_deriv_y_fill(double *where, double *where_deriv, double x, 
 		qpms_l_t lMax, gsl_sf_legendre_t lnorm, double csphase); 
 
 qpms_vswfset_sph_t *qpms_vswfset_make(qpms_l_t lMax, sph_t kdlj,
-	qpms_bessel_t btyp, qpms_normalization_t norm);//NI
+	qpms_bessel_t btyp, qpms_normalisation_t norm);//NI
-void qpms_vswfst_sph_pfree(qpms_vswfset_t *);//NI
+void qpms_vswfset_sph_pfree(qpms_vswfset_sph_t *);//NI
 
 double *qpms_legendre_y_get(double x, qpms_l_t lMax, qpms_normalisation_t norm);//NI
-double *qpms_legendre0d_y_get(qpms_l_t lMax, qpms_normalization_t norm); //NI
+double *qpms_legendre0d_y_get(qpms_l_t lMax, qpms_normalisation_t norm); //NI
-double *qpms_legendre_plus1d_y_get(qpms_l_t lMax, qpms_normalization_t norm); //NI
+double *qpms_legendre_plus1d_y_get(qpms_l_t lMax, qpms_normalisation_t norm); //NI
-double *qpms_legendre_minus1d_y_get(qpms_l_t lMax, qpms_normalization_t norm); //NI
+double *qpms_legendre_minus1d_y_get(qpms_l_t lMax, qpms_normalisation_t norm); //NI
 
 
 
-// array of pi, tau auxillary function (see [1,(37)])
+// array of Legendre and pi, tau auxillary functions (see [1,(37)])
+// This should handle correct evaluation for theta -> 0 and theta -> pi
 typedef struct {
-	//qpms_normalization_t norm;
+	//qpms_normalisation_t norm;
 	qpms_l_t lMax;
 	//qpms_y_t nelem;
-	double *pi, *tau;
+	double *leg, *pi, *tau;
 } qpms_pitau_t;
+qpms_pitau_t qpms_pitau_get(double theta, qpms_l_t lMax, qpms_normalisation_t norm);
 void qpms_pitau_free(qpms_pitau_t);//NI
 void qpms_pitau_pfree(qpms_pitau_t*);//NI