Slatec(3pm)

📦

version 2.085 · 4 versions

💖 Slatec 💖

the PDL::Slatec manual, part 1 of 3 — bestie edition 🌸

okay so real talk, this is Slatec(3pm) from the User Contributed Perl Documentation section, and it's a whole vibe. let's get into it 💅

✨ NAME ✨

PDL::Slatec — the PDL interface to the slatec numerical programming library 🧚

🌈 SYNOPSIS 🌈

use PDL::Slatec;

($ndeg, $r, $ierr, $c) = polyfit($x, $y, $w, $maxdeg, $eps);

🦋 DESCRIPTION 🦋

so this module is doing two whole jobs at once, main character behavior fr:

it gives you an interface to parts of the slatec library
it also shows you how to interface PDL to an external library, like a little tutorial bestie 🫶

using this library requires a fortran compiler — the source for the routines is provided for convenience.

right now, what's available: routines to manipulate matrices 🧮, calculate FFT's ⚡, fit data using polynomials 📈, and interpolate/integrate data using piecewise cubic Hermite interpolation ✨.

🎀 Piecewise cubic Hermite interpolation (PCHIP) 🎀

PCHIP is the slatec package of routines that do piecewise cubic Hermite interpolation of data. it comes with software to produce a monotone and "visually pleasing" interpolant to monotone data. and per Fritsch & Carlson ("Monotone piecewise cubic interpolation", SIAM Journal on Numerical Analysis 17, 2 (April 1980), pp. 238-246), that kind of interpolant can be more reasonable than a cubic spline if your data has both "steep" and "flat" sections. interpolating cumulative probability distribution functions? also a whole application, slay. these routines are cryptically named (blame FORTRAN 💀), they all start with 'ch', and they accept either float or double ndarrays.

most of the routines need an integer parameter called "check" — if you set it to 0, then no checks on the validity of the input data are made, otherwise those checks do happen. and you can set "check" to 0 if a routine like "chim" has already been successfully called 👑

so the workflow is giving:

If not known, estimate derivative values for the points using the "chim", "chic", or "chsp" routines (the following routines require both the function ("f") and derivative ("d") values at a set of points ("x")).
Evaluate, integrate, or differentiate the resulting PCH function using the routines: "chfd"; "chfe"; "chia"; "chid".
If desired, you can check the monotonicity of your data using "chcm".

👑 FUNCTIONS 👑

⭐ eigsys ⭐

Eigenvalues and eigenvectors of a real positive definite symmetric matrix.

($eigvals,$eigvecs) = eigsys($mat)

Note: this function should be extended to calculate only eigenvalues if called in scalar context! (not me calling out a TODO in the docs 😭)

⭐ matinv ⭐

Inverse of a square matrix.

($inv) = matinv($mat)

⭐ polyfit ⭐

a convenience wrapper routine hanging around the "polfit" "slatec" function. it separates supplied arguments and return values, so you don't have to deal with all that mess yourself 💅

Fit discrete data in a least squares sense by polynomials in one variable. it handles broadcasting correctly — you can literally pass in a 2D PDL (as $y) and it will pass back a 2D PDL, the rows of which are the polynomial regression results (in $r corresponding to the rows of $y.

($ndeg, $r, $ierr, $c, $coeffs, $rms) = polyfit($x, $y, $w, $maxdeg, [$eps]);

$coeffs = polyfit($x,$y,$w,$maxdeg,[$eps]);

where on input:

$x and $y are the values to fit to a polynomial. $w are weighting factors. $maxdeg is the maximum degree of polynomial to use and $eps is the required degree of fit.

and the output switches on list/scalar context 🎀

In list context:

$ndeg is the degree of polynomial actually used. $r is the values of the fitted polynomial. $ierr is a return status code, and $c is some working array or other (preserved for historical purposes). $coeffs is the polynomial coefficients of the best fit polynomial. $rms is the rms error of the fit.

In scalar context, only $coeffs is returned.

Historically, $eps was modified in-place to be a return value of the rms error. This usage is deprecated, and $eps is an optional parameter now. It is still modified if present. 🙈

$c is a working array accessible to Slatec — you can feed it to several other Slatec routines to get nice things out. It does not broadcast correctly and should probably be fixed by someone. If you are reading this, that someone might be you. 👀

This version of polyfit handles bad values correctly. Bad values in $y are ignored for the fit and give computed values on the fitted curve in the return. Bad values in $x or $w are ignored for the fit and result in bad elements in the output.

⭐ polycoef ⭐

a convenience wrapper routine around the "pcoef" "slatec" function. separates supplied arguments and return values ✨

Convert the "polyfit"/"polfit" coefficients to Taylor series form.

$tc = polycoef($l, $c, $x);

⭐ polyvalue ⭐

a convenience wrapper routine around the "pvalue" "slatec" function. separates supplied arguments and return values ✨

For multiple input x positions, a corresponding y position is calculated.

The derivatives PDL is one dimensional (of size "nder") if a single x position is supplied, two dimensional if more than one x position is supplied.

Use the coefficients "c" generated by "polyfit" (or "polfit") to evaluate the polynomial fit of degree "l", along with the first "nder" of its derivatives, at a specified point "x".

($yfit, $yp) = polyvalue($l, $nder, $x, $c);

⭐ detslatec ⭐

compute the determinant of an invertible matrix 🔥

$mat = zeroes(5,5); $mat->diagonal(0,1) .= 1; # unity matrix
$det = detslatec $mat;

Usage:

$determinant = detslatec $matrix;

Signature: detslatec(mat(n,m); [o] det())

"detslatec" computes the determinant of an invertible matrix and barfs if the matrix argument provided is non-invertible. The matrix broadcasts as usual.

This routine was previously known as "det" which clashes now with det which is provided by PDL::MatrixOps. so yeah, glow-up + rename, iconic 💅

⭐ fft ⭐

Fast Fourier Transform ⚡

$v_in = pdl(1,0,1,0);
($azero,$x,$y) = PDL::Slatec::fft($v_in);

"PDL::Slatec::fft" is a convenience wrapper for "ezfftf", and performs a Fast Fourier Transform on an input vector $v_in. The return values are the same as for "ezfftf".

⭐ rfft ⭐

reverse Fast Fourier Transform 🔁

$v_out = PDL::Slatec::rfft($azero,$x,$y);
print $v_in, $vout
[1 0 1 0] [1 0 1 0]

"PDL::Slatec::rfft" is a convenience wrapper for "ezfftb", and performs a reverse Fast Fourier Transform. The input is the same as the output of "PDL::Slatec::fft", and the output of "rfft" is a data vector, similar to what is input into "PDL::Slatec::fft".

⭐ svdc ⭐

Signature: (x(n,p);[o]s(p);[o]e(p);[o]u(n,p);[o]v(p,p);[o]work(n);longlong job();longlong [o]info())

singular value decomposition of a matrix 🧮

svdc does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ poco ⭐

Signature: (a(n,n);rcond();[o]z(n);longlong [o]info())

Factor a real symmetric positive definite matrix and estimate the condition number of the matrix.

poco does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ geco ⭐

Signature: (a(n,n);longlong [o]ipvt(n);[o]rcond();[o]z(n))

Factor a matrix using Gaussian elimination and estimate the condition number of the matrix.

geco does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ gefa ⭐

Signature: (a(n,n);longlong [o]ipvt(n);longlong [o]info())

Factor a matrix using Gaussian elimination.

gefa does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ podi ⭐

Signature: (a(n,n);[o]det(two=2);longlong job())

Compute the determinant and inverse of a certain real symmetric positive definite matrix using the factors computed by "poco".

podi does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ gedi ⭐

Signature: (a(n,n);longlong [o]ipvt(n);[o]det(two=2);[o]work(n);longlong job())

Compute the determinant and inverse of a matrix using the factors computed by "geco" or "gefa".

gedi does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ gesl ⭐

Signature: (a(lda,n);longlong ipvt(n);b(n);longlong job())

Solve the real system "A*X=B" or "TRANS(A)*X=B" using the factors computed by "geco" or "gefa".

gesl does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ rs ⭐

Signature: (a(n,n);[o]w(n);longlong matz();[o]z(n,n);[t]fvone(n);[t]fvtwo(n);longlong [o]ierr())

This subroutine calls the recommended sequence of subroutines from the eigensystem subroutine package (EISPACK) to find the eigenvalues and eigenvectors (if desired) of a REAL SYMMETRIC matrix.

rs does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ ezffti ⭐

Signature: (longlong n();[o]wsave(foo))

Subroutine ezffti initializes the work array wsave() which is used in both "ezfftf" and "ezfftb". The prime factorization of "n" together with a tabulation of the trigonometric functions are computed and stored in wsave().

ezffti does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ ezfftf ⭐

Signature: (r(n);[o]azero();[o]a(n);[o]b(n);wsave(foo))

ezfftf does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ ezfftb ⭐

Signature: ([o]r(n);azero();a(n);b(n);wsave(foo))

ezfftb does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ pcoef ⭐

Signature: (longlong l();c();[o]tc(bar);a(foo))

Convert the "polfit" coefficients to Taylor series form. "c" and a() must be of the same type.

pcoef does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ pvalue ⭐

Signature: (longlong l();x();[o]yfit();[o]yp(nder);a(foo))

Use the coefficients generated by "polfit" to evaluate the polynomial fit of degree "l", along with the first "nder" of its derivatives, at a specified point. "x" and "a" must be of the same type.

pvalue does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays.

⭐ chim ⭐

Signature: (x(n);f(n);[o]d(n);longlong [o]ierr())

Calculate the derivatives of (x,f(x)) using cubic Hermite interpolation.

Calculate the derivatives at the given set of points ("$x,$f", where $x is strictly increasing). The resulting set of points - "$x,$f,$d", referred to as the cubic Hermite representation - can then be used in other functions, such as "chfe", "chfd", and "chia".

🌸 okay that's part 1 bestie, more coming, don't go anywhere 🌸

The boundary conditions are lowkey compatible with monotonicity, bestie. And if your data is only piecewise monotonic, the interpolant will have an extremum at the switch points. For more control over these issues? Use chic — she's the upgrade. 💅

Error status returned by $ierr:

0 if successful.
> 0 if there were ierr switches in the direction of monotonicity (data still valid).
-1 if nelem($x) < 2.
-3 if $x is not strictly increasing.

And bestie, chim does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. So like, keep your data clean. 🧼

✨ chic ✨

Signature: (longlong ic(two=2);vc(two=2);mflag();x(n);f(n);[o]d(n);wk(nwk);longlong [o]ierr())

So chic calculates the derivatives of (x,f(x)) using cubic Hermite interpolation. It's giving math but make it cute. 🎀

Basically, you calculate the derivatives at the given points ($x,$f, where $x is strictly increasing). Control over the boundary conditions is given by the $ic and $vc ndarrays, and the value of $mflag determines the treatment of points where monotonicity switches direction. A simpler, more restricted, interface is available using chim.

The first and second elements of $ic determine the boundary conditions at the start and end of the data respectively. If the value is 0, then the default condition, as used by chim, is adopted. If greater than zero, no adjustment for monotonicity is made, otherwise if less than zero the derivative will be adjusted. The allowed magnitudes for ic(0) are:

1 if first derivative at x(0) is given in vc(0).
2 if second derivative at x(0) is given in vc(0).
3 to use the 3-point difference formula for d(0). (Reverts to the default b.c. if n < 3)
4 to use the 4-point difference formula for d(0). (Reverts to the default b.c. if n < 4)
5 to set d(0) so that the second derivative is continuous at x(1). (Reverts to the default b.c. if n < 4)

The values for ic(1) are the same as above, except that the first-derivative value is stored in vc(1) for cases 1 and 2. The values of $vc need only be set if options 1 or 2 are chosen for $ic.

Set $mflag = 0 if interpolant is required to be monotonic in each interval, regardless of the data. This causes $d to be set to 0 at all switch points. Set $mflag to be non-zero to use a formula based on the 3-point difference formula at switch points. If $mflag > 0, then the interpolant at switch points is forced to not deviate from the data by more than $mflag*dfloc, where dfloc is the maximum of the change of $f on this interval and its two immediate neighbours. If $mflag < 0, no such control is to be imposed.

The ndarray $wk is only needed for work space. However, I could not get it to work as a temporary variable, so you must supply it; it is a 1D ndarray with 2*n elements. Like, it's a whole thing. 🙈

Error status returned by $ierr:

0 if successful.
1 if ic(0) < 0 and d(0) had to be adjusted for monotonicity.
2 if ic(1) < 0 and d(n-1) had to be adjusted for monotonicity.
3 if both 1 and 2 are true.
-1 if n < 2.
-3 if $x is not strictly increasing.
-4 if abs(ic(0)) > 5.
-5 if abs(ic(1)) > 5.
-6 if both -4 and -5 are true.
-7 if nwk < 2*(n-1).

And bestie, chic does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chsp ✨

Signature: (longlong ic(two=2);vc(two=2);x(n);f(n);[o]d(n);wk(nwk);longlong [o]ierr())

So chsp calculates the derivatives of (x,f(x)) using cubic spline interpolation. It's like chic but with splines, slay. 💅

Calculate the derivatives, using cubic spline interpolation, at the given points ($x,$f), with the specified boundary conditions. Control over the boundary conditions is given by the $ic and $vc ndarrays. The resulting values — $x,$f,$d — can be used in all the functions which expect a cubic Hermite function.

The first and second elements of $ic determine the boundary conditions at the start and end of the data respectively. The allowed values for ic(0) are:

0 to set d(0) so that the third derivative is continuous at x(1).
1 if first derivative at x(0) is given in vc(0).
2 if second derivative at "x(0") is given in vc(0).
3 to use the 3-point difference formula for d(0). (Reverts to the default b.c. if n < 3.)
4 to use the 4-point difference formula for d(0). (Reverts to the default b.c. if n < 4.)

The ndarray $wk is only needed for work space. However, I could not get it to work as a temporary variable, so you must supply it; it is a 1D ndarray with 2*n elements. Same energy as chic, bestie. 🙈

Error status returned by $ierr:

0 if successful.
-1 if nelem($x) < 2.
-3 if $x is not strictly increasing.
-4 if ic(0) < 0 or ic(0) > 4.
-5 if ic(1) < 0 or ic(1) > 4.
-6 if both of the above are true.
-7 if nwk < 2*n.
-8 in case of trouble solving the linear system for the interior derivative values.

And chsp does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chfd ✨

Signature: (x(n);f(n);d(n);longlong check();xe(ne);[o]fe(ne);[o]de(ne);longlong [o]ierr())

Interpolate function and derivative values. It's giving double duty. 👑

Given a piecewise cubic Hermite function — such as from chim — evaluate the function ($fe) and derivative ($de) at a set of points ($xe). If function values alone are required, use chfe. Set check to 0 to skip checks on the input data.

Error status returned by $ierr:

0 if successful.
>0 if extrapolation was performed at ierr points (data still valid).
-1 if nelem($x) < 2
-3 if $x is not strictly increasing.
-4 if nelem($xe) < 1.
-5 if an error has occurred in a lower-level routine, which should never happen.

And chfd does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chfe ✨

Signature: (x(n);f(n);d(n);longlong check();xe(ne);[o]fe(ne);longlong [o]ierr())

Interpolate function values. Just the function, no derivative, we keep it simple. ✨

Given a piecewise cubic Hermite function — such as from chim — evaluate the function ($fe) at a set of points ($xe). If derivative values are also required, use chfd. Set check to 0 to skip checks on the input data.

Error status returned by $ierr:

0 if successful.
>0 if extrapolation was performed at ierr points (data still valid).
-1 if nelem($x) < 2
-3 if $x is not strictly increasing.
-4 if nelem($xe) < 1.

And chfe does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chia ✨

Signature: (x(n);f(n);d(n);longlong check();la();lb();[o]ans();longlong [o]ierr())

Integrate (x,f(x)) over arbitrary limits. Math but make it integrals. 📐

Evaluate the definite integral of a piecewise cubic Hermite function over an arbitrary interval, given by [$la,$lb]. $d should contain the derivative values, computed by chim. See chid if the integration limits are data points. Set check to 0 to skip checks on the input data.

The values of $la and $lb do not have to lie within $x, although the resulting integral value will be highly suspect if they are not. So like, proceed with caution, bestie. ⚠️

Error status returned by $ierr:

0 if successful.
1 if $la lies outside $x.
2 if $lb lies outside $x.
3 if both 1 and 2 are true.
-1 if nelem($x) < 2
-3 if $x is not strictly increasing.
-4 if an error has occurred in a lower-level routine, which should never happen.

And chia does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chid ✨

Signature: (x(n);f(n);d(n);longlong check();longlong ia();longlong ib();[o]ans();longlong [o]ierr())

Integrate (x,f(x)) between data points. It's like chia but with data points, iykyk. 💅

Evaluate the definite integral of a a piecewise cubic Hermite function between x($ia) and x($ib). (Yes, the original says "a a", we're keeping it, it's iconic. 😭)

See chia for integration between arbitrary limits.

Although using a fortran routine, the values of $ia and $ib are zero offset. $d should contain the derivative values, computed by chim. Set check to 0 to skip checks on the input data.

Error status returned by $ierr:

0 if successful.
-1 if nelem($x) < 2.
-3 if $x is not strictly increasing.
-4 if $ia or $ib is out of range.

And chid does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chcm ✨

Signature: (x(n);f(n);d(n);longlong check();longlong [o]ismon(n);longlong [o]ierr())

Check the given piecewise cubic Hermite function for monotonicity. Vibe check! ✨

The outout ndarray $ismon indicates over which intervals the function is monotonic. Set check to 0 to skip checks on the input data.

For the data interval [x(i),x(i+1)], the values of ismon(i) can be:

-3 if function is probably decreasing
-1 if function is strictly decreasing
0 if function is constant
1 if function is strictly increasing
2 if function is non-monotonic
3 if function is probably increasing

If abs(ismon(i)) == 3, the derivative values are near the boundary of the monotonicity region. A small increase produces non-monotonicity, whereas a decrease produces strict monotonicity. So like, tread carefully. 🚨

The above applies to i = 0 .. nelem($x)-1. The last element of $ismon indicates whether the entire function is monotonic over $x.

Error status returned by $ierr:

0 if successful.
-1 if n < 2.
-3 if $x is not strictly increasing.

And chcm does NOT process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. 🧼

✨ chbs ✨

Signature: (x(n);f(n);d(n);longlong knotyp();longlong nknots();t(tsize);[o]bcoef(bsize);longlong [o]ndim();longlong [o]kord();longlong [o]ierr())

Piecewise Cubic Hermite function to B-Spline converter. It's giving transformation. 🦋

The resulting B-spline representation of the data (i.e. nknots, t, bcoeff, ndim, and kord) can be evaluated by bvalu (which is currently not available). Sad, but okay. 😭

Array sizes: tsize = 2*n + 4, bsize = 2*n, and ndim = 2*n.

knotyp is a flag which controls the knot sequence. The knot sequence t is normally computed from $x by putting a double knot at each x and setting the end knot pairs according to the value of knotyp (where m = ndim = 2*n):

0 — quadruple knots at the first point and the last point. Slay, boundaries secured. ✨
1 — replicate the lengths of the extreme subintervals: "t( 0 ) = t( 1 ) = x(0) - (x(1)-x(0))" and "t(m+3) = t(m+2) = x(n-1) + (x(n-1)-x(n-2))" — bestie it's literally copy-paste energy from the end intervals.
2 — periodic placement of the boundary knots: "t( 0 ) = t( 1 ) = x(0) - (x(n-1)-x(n-2))" and "t(m+3) = t(m+2) = x(n) + (x(1)-x(0))" — it's giving "wraparound era". 🌀
<0 — assume the nknots and t were set in a previous call. iykyk, we're not redoing that work. 💅

nknots is the number of knots, and it may be changed by the routine. If knotyp >= 0, then nknots gets set to ndim+4. Otherwise it's an input variable, and bestie — an error will occur if its value is not equal to ndim+4. So don't play. 😤

t is the array of 2*n+4 knots for the B-representation, and it may be changed by the routine too. If knotyp >= 0, t will be changed so the interior double knots equal the x-values and the boundary knots are set as indicated above. Otherwise it's assumed that t was set by a previous call — and no cap, no check is made to verify that the data forms a legitimate knot sequence. So if you hand it nonsense, that's on you. 🫣

Error status returned by $ierr:

0 if successful. 🎉
-4 if knotyp > 2.
-5 if knotyp < 0 and nknots != 2*n + 4.

🚨 chbs does not process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. That's just how it is, fr fr.

🧮 polfit 🧮

Signature:

Signature: (x(n); y(n); w(n); longlong maxdeg(); longlong [o]ndeg(); [o]eps(); [o]r(n); longlong [o]ierr(); [o]a(foo); [o]coeffs(bar);[t]xtmp(n);[t]ytmp(n);[t]wtmp(n);[t]rtmp(n))

Fit discrete data in a least squares sense by polynomials in one variable. 🦋 And listen: x(), y() and w() must be of the same type — they're a matching set, don't break up the trio. This version handles bad values appropriately, and honestly? Obsessed. 💖

✨ polfit does process bad values. It will set the bad-value flag of all output ndarrays if the flag is set for any of the input ndarrays. The opposite of chbs energy, and we love the range.

✍️ AUTHOR ✍️

There is no warranty. 😭 You are allowed to redistribute this software / documentation under certain conditions. For details, see the file COPYING in the PDL distribution. If this file is separated from the PDL distribution, the copyright notice should be included in the file. Periodt. 💅

and that's the whole page, bestie 🫶 go forth and spline responsibly — this was Slatec(3pm) from perl v5.38.2, dated 2024-04-10. bye luv u 💜🎀

📦 where this page comes from

ubuntu noble: pdl 1:2.085-1ubuntu1 · Artistic-1.0-Perl OR GPL-1.0-or-later
💖 girlified by deepseek-ai/DeepSeek-V4.1-Flash on 2026-10-14: a language model's rewrite of the original, under its terms