Lines Matching +full:- +full:m
8 REQUIRES configuration of PETSc with option --download-adolc.
10 For documentation on ADOL-C, see
11 $PETSC_ARCH/externalpackages/ADOL-C-2.6.0/ADOL-C/doc/adolc-manual.pdf
14 /* --------------------------------------------------------------------------------
16 ----------------------------------------------------------------------------- */
21 assembled (not recommended for non-toy problems!).
24 tag - tape identifier
25 u_vec - vector at which to evaluate Jacobian
26 ctx - ADOL-C context, as defined above
29 A - Mat object corresponding to Jacobian
34 PetscInt i, j, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeRHSJacobian() local
38 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeRHSJacobian()
39 if (adctx->Seed) fov_forward(tag, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeRHSJacobian()
40 else jacobian(tag, m, n, u_vec, J); in PetscAdolcComputeRHSJacobian()
41 if (adctx->sparse) { in PetscAdolcComputeRHSJacobian()
42 PetscCall(RecoverJacobian(A, INSERT_VALUES, m, p, adctx->Rec, J, NULL)); in PetscAdolcComputeRHSJacobian()
44 for (i = 0; i < m; i++) { in PetscAdolcComputeRHSJacobian()
46 … if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValues(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeRHSJacobian()
59 assembled (not recommended for non-toy problems!).
62 tag - tape identifier
63 u_vec - vector at which to evaluate Jacobian
64 ctx - ADOL-C context, as defined above
67 A - Mat object corresponding to Jacobian
72 PetscInt i, j, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeRHSJacobianLocal() local
76 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeRHSJacobianLocal()
77 if (adctx->Seed) fov_forward(tag, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeRHSJacobianLocal()
78 else jacobian(tag, m, n, u_vec, J); in PetscAdolcComputeRHSJacobianLocal()
79 if (adctx->sparse) { in PetscAdolcComputeRHSJacobianLocal()
80 PetscCall(RecoverJacobianLocal(A, INSERT_VALUES, m, p, adctx->Rec, J, NULL)); in PetscAdolcComputeRHSJacobianLocal()
82 for (i = 0; i < m; i++) { in PetscAdolcComputeRHSJacobianLocal()
84 …if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValuesLocal(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeRHSJacobianLocal()
97 assembled (not recommended for non-toy problems!).
100 tag1 - tape identifier for df/dx part
101 tag2 - tape identifier for df/d(xdot) part
102 u_vec - vector at which to evaluate Jacobian
103 ctx - ADOL-C context, as defined above
106 A - Mat object corresponding to Jacobian
111 PetscInt i, j, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeIJacobian() local
115 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeIJacobian()
118 if (adctx->Seed) fov_forward(tag1, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobian()
119 else jacobian(tag1, m, n, u_vec, J); in PetscAdolcComputeIJacobian()
121 if (adctx->sparse) { in PetscAdolcComputeIJacobian()
122 PetscCall(RecoverJacobian(A, INSERT_VALUES, m, p, adctx->Rec, J, NULL)); in PetscAdolcComputeIJacobian()
124 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobian()
126 … if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValues(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeIJacobian()
134 if (adctx->Seed) fov_forward(tag2, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobian()
135 else jacobian(tag2, m, n, u_vec, J); in PetscAdolcComputeIJacobian()
136 if (adctx->sparse) { in PetscAdolcComputeIJacobian()
137 PetscCall(RecoverJacobian(A, ADD_VALUES, m, p, adctx->Rec, J, &a)); in PetscAdolcComputeIJacobian()
139 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobian()
141 if (fabs(J[i][j]) > 1.e-16) { in PetscAdolcComputeIJacobian()
161 tag - tape identifier for df/dx part
162 u_vec - vector at which to evaluate Jacobian
163 ctx - ADOL-C context, as defined above
166 A - Mat object corresponding to Jacobian
171 PetscInt i, j, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeIJacobianIDMass() local
175 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeIJacobianIDMass()
178 if (adctx->Seed) fov_forward(tag, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobianIDMass()
179 else jacobian(tag, m, n, u_vec, J); in PetscAdolcComputeIJacobianIDMass()
181 if (adctx->sparse) { in PetscAdolcComputeIJacobianIDMass()
182 PetscCall(RecoverJacobian(A, INSERT_VALUES, m, p, adctx->Rec, J, NULL)); in PetscAdolcComputeIJacobianIDMass()
184 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobianIDMass()
186 … if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValues(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeIJacobianIDMass()
202 assembled (not recommended for non-toy problems!).
205 tag1 - tape identifier for df/dx part
206 tag2 - tape identifier for df/d(xdot) part
207 u_vec - vector at which to evaluate Jacobian
208 ctx - ADOL-C context, as defined above
211 A - Mat object corresponding to Jacobian
216 PetscInt i, j, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeIJacobianLocal() local
220 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeIJacobianLocal()
223 if (adctx->Seed) fov_forward(tag1, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobianLocal()
224 else jacobian(tag1, m, n, u_vec, J); in PetscAdolcComputeIJacobianLocal()
225 if (adctx->sparse) { in PetscAdolcComputeIJacobianLocal()
226 PetscCall(RecoverJacobianLocal(A, INSERT_VALUES, m, p, adctx->Rec, J, NULL)); in PetscAdolcComputeIJacobianLocal()
228 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobianLocal()
230 …if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValuesLocal(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeIJacobianLocal()
238 if (adctx->Seed) fov_forward(tag2, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobianLocal()
239 else jacobian(tag2, m, n, u_vec, J); in PetscAdolcComputeIJacobianLocal()
240 if (adctx->sparse) { in PetscAdolcComputeIJacobianLocal()
241 PetscCall(RecoverJacobianLocal(A, ADD_VALUES, m, p, adctx->Rec, J, &a)); in PetscAdolcComputeIJacobianLocal()
243 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobianLocal()
245 if (fabs(J[i][j]) > 1.e-16) { in PetscAdolcComputeIJacobianLocal()
265 tag - tape identifier for df/dx part
266 u_vec - vector at which to evaluate Jacobian
267 ctx - ADOL-C context, as defined above
270 A - Mat object corresponding to Jacobian
275 PetscInt i, j, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeIJacobianLocalIDMass() local
279 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeIJacobianLocalIDMass()
282 if (adctx->Seed) fov_forward(tag, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobianLocalIDMass()
283 else jacobian(tag, m, n, u_vec, J); in PetscAdolcComputeIJacobianLocalIDMass()
284 if (adctx->sparse) { in PetscAdolcComputeIJacobianLocalIDMass()
285 PetscCall(RecoverJacobianLocal(A, INSERT_VALUES, m, p, adctx->Rec, J, NULL)); in PetscAdolcComputeIJacobianLocalIDMass()
287 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobianLocalIDMass()
289 …if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValuesLocal(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeIJacobianLocalIDMass()
302 /* --------------------------------------------------------------------------------
304 ----------------------------------------------------------------------------- */
310 tag - tape identifier
311 u_vec - vector at which to evaluate Jacobian
312 params - the parameters w.r.t. which we differentiate
313 ctx - ADOL-C context, as defined above
316 A - Mat object corresponding to Jacobian
321 PetscInt i, j = 0, m = adctx->m, n = adctx->n, p = adctx->num_params; in PetscAdolcComputeRHSJacobianP() local
326 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeRHSJacobianP()
334 fov_forward(tag, m, n + p, p, concat, S, NULL, J); in PetscAdolcComputeRHSJacobianP()
339 for (i = 0; i < m; i++) { in PetscAdolcComputeRHSJacobianP()
341 if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValues(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeRHSJacobianP()
354 tag - tape identifier
355 u_vec - vector at which to evaluate Jacobian
356 params - the parameters w.r.t. which we differentiate
357 ctx - ADOL-C context, as defined above
360 A - Mat object corresponding to Jacobian
365 PetscInt i, j = 0, m = adctx->m, n = adctx->n, p = adctx->num_params; in PetscAdolcComputeRHSJacobianPLocal() local
370 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeRHSJacobianPLocal()
378 fov_forward(tag, m, n + p, p, concat, S, NULL, J); in PetscAdolcComputeRHSJacobianPLocal()
383 for (i = 0; i < m; i++) { in PetscAdolcComputeRHSJacobianPLocal()
385 …if (fabs(J[i][j]) > 1.e-16) PetscCall(MatSetValuesLocal(A, 1, &i, 1, &j, &J[i][j], INSERT_VALUES)); in PetscAdolcComputeRHSJacobianPLocal()
394 /* --------------------------------------------------------------------------------
396 ----------------------------------------------------------------------------- */
403 tag1 - tape identifier for df/dx part
404 tag2 - tape identifier for df/d(xdot) part
405 u_vec - vector at which to evaluate Jacobian
406 ctx - ADOL-C context, as defined above
409 diag - Vec object corresponding to Jacobian diagonal
414 PetscInt i, m = adctx->m, n = adctx->n, p = adctx->p; in PetscAdolcComputeIJacobianAndDiagonalLocal() local
418 PetscCall(AdolcMalloc2(m, p, &J)); in PetscAdolcComputeIJacobianAndDiagonalLocal()
421 if (adctx->Seed) fov_forward(tag1, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobianAndDiagonalLocal()
422 else jacobian(tag1, m, n, u_vec, J); in PetscAdolcComputeIJacobianAndDiagonalLocal()
423 if (adctx->sparse) { in PetscAdolcComputeIJacobianAndDiagonalLocal()
424 PetscCall(RecoverDiagonalLocal(diag, INSERT_VALUES, m, adctx->rec, J, NULL)); in PetscAdolcComputeIJacobianAndDiagonalLocal()
426 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobianAndDiagonalLocal()
427 … if (fabs(J[i][i]) > 1.e-16) PetscCall(VecSetValuesLocal(diag, 1, &i, &J[i][i], INSERT_VALUES)); in PetscAdolcComputeIJacobianAndDiagonalLocal()
434 if (adctx->Seed) fov_forward(tag2, m, n, p, u_vec, adctx->Seed, NULL, J); in PetscAdolcComputeIJacobianAndDiagonalLocal()
435 else jacobian(tag2, m, n, u_vec, J); in PetscAdolcComputeIJacobianAndDiagonalLocal()
436 if (adctx->sparse) { in PetscAdolcComputeIJacobianAndDiagonalLocal()
437 PetscCall(RecoverDiagonalLocal(diag, ADD_VALUES, m, adctx->rec, J, NULL)); in PetscAdolcComputeIJacobianAndDiagonalLocal()
439 for (i = 0; i < m; i++) { in PetscAdolcComputeIJacobianAndDiagonalLocal()
440 if (fabs(J[i][i]) > 1.e-16) { in PetscAdolcComputeIJacobianAndDiagonalLocal()