xref: /libCEED/tests/t301-basis.c (revision 9ac7b42ea0005d14dd8105f72e4e8b70779a82a1)
14411cf47Sjeremylt /// @file
252bfb9bbSJeremy L Thompson /// Test QR Factorization
352bfb9bbSJeremy L Thompson /// \test Test QR Factorization
457c64913Sjeremylt #include <ceed.h>
5*9ac7b42eSJeremy L Thompson #include <ceed/backend.h>
6*9ac7b42eSJeremy L Thompson #include <math.h>
757c64913Sjeremylt 
857c64913Sjeremylt int main(int argc, char **argv) {
957c64913Sjeremylt   Ceed ceed;
10*9ac7b42eSJeremy L Thompson   CeedScalar A[12]    = {1, -1, 4, 1, 4, -2, 1, 4, 2, 1, -1, 0};
1152bfb9bbSJeremy L Thompson   CeedScalar qr[12]   = {1, -1, 4, 1, 4, -2, 1, 4, 2, 1, -1, 0};
12*9ac7b42eSJeremy L Thompson   CeedScalar A_qr[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
13*9ac7b42eSJeremy L Thompson   CeedScalar tau[4];
1457c64913Sjeremylt 
1557c64913Sjeremylt   CeedInit(argv[1], &ceed);
16aedaa0e5Sjeremylt 
1752bfb9bbSJeremy L Thompson   CeedQRFactorization(ceed, qr, tau, 4, 3);
18*9ac7b42eSJeremy L Thompson   for (CeedInt i=0; i<3; i++)
19*9ac7b42eSJeremy L Thompson     for (CeedInt j=i; j<3; j++)
20*9ac7b42eSJeremy L Thompson       A_qr[i*3+j] = qr[i*3+j];
21*9ac7b42eSJeremy L Thompson   CeedHouseholderApplyQ(A_qr, qr, tau, CEED_NOTRANSPOSE, 4, 3, 3, 3, 1);
22*9ac7b42eSJeremy L Thompson 
23*9ac7b42eSJeremy L Thompson   for (CeedInt i=0; i<12; i++)
24*9ac7b42eSJeremy L Thompson     if (fabs(A_qr[i] - A[i]) > 10*CEED_EPSILON)
25*9ac7b42eSJeremy L Thompson       // LCOV_EXCL_START
26*9ac7b42eSJeremy L Thompson       printf("Error in QR factorization A_qr[%d] = %f != A[%d] = %f\n",
27*9ac7b42eSJeremy L Thompson              i, A_qr[i], i, A[i]);
28*9ac7b42eSJeremy L Thompson   // LCOV_EXCL_STOP
29*9ac7b42eSJeremy L Thompson 
3057c64913Sjeremylt   CeedDestroy(&ceed);
3157c64913Sjeremylt   return 0;
3257c64913Sjeremylt }
33