Skip to content

Instantly share code, notes, and snippets.

@staticfloat
Created April 5, 2023 15:12
Show Gist options
  • Save staticfloat/494272e86b44a3e2d31de898473b3c64 to your computer and use it in GitHub Desktop.
Save staticfloat/494272e86b44a3e2d31de898473b3c64 to your computer and use it in GitHub Desktop.
Accelerate NEWLAPCK dpstrf bug
#include <stdio.h>
#include <stdint.h>
#ifdef ILP64
typedef int64_t blasint;
#else
typedef int32_t blasint;
#endif
extern blasint dpstrf(char *, blasint *, double *, blasint *, blasint *, blasint *, double *, double*, long*);
#define N 4
int main()
{
blasint pivot[N];
long info;
double A[N][N];
double work[2*N];
blasint order = N;
blasint rank = 0;
blasint lda = N;
blasint stride = 1;
double tol = 0.0;
// Initialize `A` with known values (transposed into FORTRAN ordering)
A[0][0] = 3.4291134; A[1][0] = -0.61112815; A[2][0] = 0.8155207; A[3][0] = -0.9183448;
A[0][1] = -0.61112815; A[1][1] = 3.1250076; A[2][1] = -0.44400603; A[3][1] = 0.97749346;
A[0][2] = 0.8155207; A[1][2] = -0.44400603; A[2][2] = 0.5413656; A[3][2] = 0.53000593;
A[0][3] = -0.9183448; A[1][3] = 0.97749346; A[2][3] = 0.53000593; A[3][3] = 5.108155;
// find solution using LAPACK routine dpstrf, all the arguments have to
// be pointers and you have to add an underscore to the routine name
dpstrf("U", &order, &A[0][0], &lda, &pivot[0], &rank, &tol, &work[0], &info);
if (info != 0) {
printf("ERROR: info == %ld!\n", info);
return 1;
}
// Print out diagonal of A
printf("diag(A):");
for (blasint j=0; j<N; j++) {
printf(" %8.4f", A[j][j]);
}
printf("\n");
return 0;
}
all: run
dpstrf_newlapack: dpstrf_test.c
$(CC) -o $@ -Ddpstrf='dpstrf$$NEWLAPACK' dpstrf_test.c -framework Accelerate
dpstrf_newlapack_ilp64: dpstrf_test.c
$(CC) -o $@ -Ddpstrf='dpstrf$$NEWLAPACK$$ILP64' -DILP64 dpstrf_test.c -framework Accelerate
dpstrf_oldlapack: dpstrf_test.c
$(CC) -o $@ -Ddpstrf='dpstrf_' dpstrf_test.c -framework Accelerate
run: dpstrf_newlapack dpstrf_newlapack_ilp64 dpstrf_oldlapack
./dpstrf_newlapack
./dpstrf_newlapack_ilp64
./dpstrf_oldlapack
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment