Slatec functions almost always will require that you allocate them a linear
array of floating points and/or integers, named something along the
lines of work[] which are needed to use as scratch space for its work. The size of this array depends on the
Slatec function and on the specific problem. Each Slatec function has the
relation for its size given in the documentation. For instance, in the example
in the previous section, cpqr79 had the relation,
C WORK REAL work array of dimension at least 2*NDEG*(NDEG+1)
You can simply give the Slatec function an array of the above size or larger. For quick and dirty programming it is not worthwhile nor neccessary to give much thought to work sizes--just give a generous, rough estimate of the size; even if you estimate turns out on the low side, Slatec is very good about printing out an error message asking for a larger work array.
If more tedious programming is for some reason required, and if you do not
want to deal with these work relations all the time,
you might want to adapt the Slatec function once and for all. This can be done
quite generically for all Slatec functions; below we demonstrate this on the
function cpqr79.
/* file: cpqr.c
cpqr(..) calls the slatec fortran function cpqr79(..)
cpqr(..) now automatically deals with memory allocation.
cpqr(..) now is also call by value for all arguments
except for the array of coefficients and roots.
The error flag is no longer an argument.
The error flag is returned when cpqr(..) exits.
***
Formula for min dimension of work[] --ref. cpqrt79.f code comment:
C WORK REAL work array of dimension at least 2*NDEG*(NDEG+1)
*/
#include<stdlib.h>
#define REAL float
struct complex {REAL re; REAL im;};
void cpqr79(int*, struct complex*, struct complex*,int*,REAL*);
REAL CPQR_STD_WK[20200]; /* Sets up a default matrix. NDEG<=100 */
int cpqr(int ndeg, struct complex* coeff, struct complex* root)
{
int error_code=0;
unsigned int workbytes;
REAL* work;
work=CPQR_STD_WK;
workbytes=sizeof(REAL)*2*ndeg*(ndeg+1);
if (workbytes>sizeof(REAL)*20200) work= (REAL*) malloc(workbytes);
cpqr79(&ndeg, coeff, root, &error_code, work);
if (workbytes>sizeof(REAL)*20200) free(work);
return(error_code); /* returns error flag */
}
We can now much more easily use the function cpqr79 by calling its
Compare the example
below with the listing in section 2.:
/* file: use.c */
#include <stdio.h>
#include <math.h>
#define REAL float
struct complex {REAL re; REAL im;};
void cpqr(int, struct complex*, struct complex*);
int main() {
int i,ndeg=2;
struct complex coeff[3];
struct complex root[2];
REAL f;
coeff[0].re=1.; coeff[0].im=0.;
coeff[1].re=3.; coeff[1].im=0.;
coeff[2].re=2.; coeff[2].im=0.;
cpqr(ndeg, coeff, root);
for (i=0;i<ndeg;i++ ) printf("(%g,%g).",root[i].re,root[i].im);
puts("\n"); return 0;
}