4 #include "t541-operator.h"
59   // QFunction - setup diff  in main()
65   // Operator - setup diff  in main()
74   // QFunction - apply  in main()
80   // Operator - apply  in main()
90   CeedScalar S[16];  in main()  local
101         u_array[p - 1] = 1.0;  in main()
104         u_array[p * p - p] = 1.0;  in main()
107         u_array[p * p - 1] = 1.0;  in main()
117     S[0 * 4 + i] = -v_array[0];  in main()
118     S[1 * 4 + i] = -v_array[p - 1];  in main()
119     S[2 * 4 + i] = -v_array[p * p - p];  in main()
120     S[3 * 4 + i] = -v_array[p * p - 1];  in main()
126 …S_inv[0] = S[5] * S[10] * S[15] - S[5] * S[11] * S[14] - S[9] * S[6] * S[15] + S[9] * S[7] * S[14]…  in main()
128 …S_inv[4] = -S[4] * S[10] * S[15] + S[4] * S[11] * S[14] + S[8] * S[6] * S[15] - S[8] * S[7] * S[14…  in main()
130 …S_inv[8] = S[4] * S[9] * S[15] - S[4] * S[11] * S[13] - S[8] * S[5] * S[15] + S[8] * S[7] * S[13] …  in main()
132 …S_inv[12] = -S[4] * S[9] * S[14] + S[4] * S[10] * S[13] + S[8] * S[5] * S[14] - S[8] * S[6] * S[13…  in main()
134 …S_inv[1] = -S[1] * S[10] * S[15] + S[1] * S[11] * S[14] + S[9] * S[2] * S[15] - S[9] * S[3] * S[14…  in main()
136 …S_inv[5] = S[0] * S[10] * S[15] - S[0] * S[11] * S[14] - S[8] * S[2] * S[15] + S[8] * S[3] * S[14]…  in main()
138 …S_inv[9] = -S[0] * S[9] * S[15] + S[0] * S[11] * S[13] + S[8] * S[1] * S[15] - S[8] * S[3] * S[13]…  in main()
140 …S_inv[13] = S[0] * S[9] * S[14] - S[0] * S[10] * S[13] - S[8] * S[1] * S[14] + S[8] * S[2] * S[13]…  in main()
142 …S_inv[2] = S[1] * S[6] * S[15] - S[1] * S[7] * S[14] - S[5] * S[2] * S[15] + S[5] * S[3] * S[14] +…  in main()
144 …S_inv[6] = -S[0] * S[6] * S[15] + S[0] * S[7] * S[14] + S[4] * S[2] * S[15] - S[4] * S[3] * S[14] …  in main()
146 …S_inv[10] = S[0] * S[5] * S[15] - S[0] * S[7] * S[13] - S[4] * S[1] * S[15] + S[4] * S[3] * S[13] …  in main()
148 …S_inv[14] = -S[0] * S[5] * S[14] + S[0] * S[6] * S[13] + S[4] * S[1] * S[14] - S[4] * S[2] * S[13]…  in main()
150 …S_inv[3] = -S[1] * S[6] * S[11] + S[1] * S[7] * S[10] + S[5] * S[2] * S[11] - S[5] * S[3] * S[10] …  in main()
152 …S_inv[7] = S[0] * S[6] * S[11] - S[0] * S[7] * S[10] - S[4] * S[2] * S[11] + S[4] * S[3] * S[10] +…  in main()
154 …S_inv[11] = -S[0] * S[5] * S[11] + S[0] * S[7] * S[9] + S[4] * S[1] * S[11] - S[4] * S[3] * S[9] -…  in main()
156 …S_inv[15] = S[0] * S[5] * S[10] - S[0] * S[6] * S[9] - S[4] * S[1] * S[10] + S[4] * S[2] * S[9] + …  in main()
158     det = 1 / (S[0] * S_inv[0] + S[1] * S_inv[4] + S[2] * S_inv[8] + S[3] * S_inv[12]);  in main()
171 …r (CeedInt j = 0; j < p; j++) u_array[i * p + j] = -(nodes[i] - 1.0) * (nodes[i] + 1.0) - (nodes[j…  in main()
181     // -- Zero corners  in main()
186     v_array[p - 1]     = 0.0;  in main()
187     v_array[p * p - p] = 0.0;  in main()
188     v_array[p * p - 1] = 0.0;  in main()
191     // -- Apply FDM inverse to interior  in main()
194     // -- Pick off corners  in main()
200     w_Pi[1] = w_array[p - 1];  in main()
201     w_Pi[2] = w_array[p * p - p];  in main()
202     w_Pi[3] = w_array[p * p - 1];  in main()
205     // -- Apply inverse of Schur complement  in main()
215     // -- Set corners  in main()
218     v_array[p - 1]     = v_Pi[1];  in main()
219     v_array[p * p - p] = v_Pi[2];  in main()
220     v_array[p * p - 1] = v_Pi[3];  in main()
223     // -- Apply full FDM inverse again  in main()
234         if (fabs(u_array[i * p + j] - w_array[i * p + j]) > 2e-3) {  in main()