Cantera  3.1.0a1
ctlapack.h
Go to the documentation of this file.
1 /**
2  * @file ctlapack.h
3  */
4 
5 // This file is part of Cantera. See License.txt in the top-level directory or
6 // at https://cantera.org/license.txt for license and copyright information.
7 
8 #ifndef CT_CTLAPACK_H
9 #define CT_CTLAPACK_H
10 
11 #include "cantera/base/ct_defs.h"
12 
13 // map BLAS names to names with or without a trailing underscore.
14 #ifndef LAPACK_FTN_TRAILING_UNDERSCORE
15 
16 #define _DGEMV_ dgemv
17 #define _DGETRF_ dgetrf
18 #define _DGETRS_ dgetrs
19 #define _DGETRI_ dgetri
20 #define _DGELSS_ dgelss
21 #define _DGBCON_ dgbcon
22 #define _DGBSV_ dgbsv
23 #define _DGBTRF_ dgbtrf
24 #define _DGBTRS_ dgbtrs
25 #define _DGECON_ dgecon
26 #define _DLANGE_ dlange
27 #define _DSCAL_ dscal
28 #define _DGEQRF_ dgeqrf
29 #define _DORMQR_ dormqr
30 #define _DTRTRS_ dtrtrs
31 #define _DTRCON_ dtrcon
32 #define _DPOTRF_ dpotrf
33 #define _DPOTRS_ dpotrs
34 
35 #else
36 
37 #define _DGEMV_ dgemv_
38 #define _DGETRF_ dgetrf_
39 #define _DGETRS_ dgetrs_
40 #define _DGETRI_ dgetri_
41 #define _DGELSS_ dgelss_
42 #define _DGBCON_ dgbcon_
43 #define _DGBSV_ dgbsv_
44 #define _DGBTRF_ dgbtrf_
45 #define _DGBTRS_ dgbtrs_
46 #define _DGECON_ dgecon_
47 #define _DLANGE_ dlange_
48 #define _DSCAL_ dscal_
49 #define _DGEQRF_ dgeqrf_
50 #define _DORMQR_ dormqr_
51 #define _DTRTRS_ dtrtrs_
52 #define _DTRCON_ dtrcon_
53 #define _DPOTRF_ dpotrf_
54 #define _DPOTRS_ dpotrs_
55 
56 #endif
57 
58 
59 namespace ctlapack
60 {
61 typedef enum {Transpose = 1, NoTranspose = 0} transpose_t;
62 typedef enum {ColMajor = 1, RowMajor = 0} storage_t;
63 typedef enum {UpperTriangular = 0, LowerTriangular = 1} upperlower_t;
64 typedef enum {Left = 0, Right = 1} side_t;
65 }
66 const char no_yes[2] = {'N', 'T'};
67 const char upper_lower[2] = {'U', 'L'};
68 const char left_right[2] = {'L', 'R'};
69 
70 // C interfaces for Fortran Lapack routines
71 extern "C" {
72 
73 #ifdef LAPACK_FTN_STRING_LEN_AT_END
74  int _DGEMV_(const char* transpose,
75  const integer* m, const integer* n, const doublereal* alpha,
76  const doublereal* a, const integer* lda, const doublereal* x,
77  const integer* incX, const doublereal* beta, doublereal* y,
78  const integer* incY, ftnlen trsize);
79 #else
80  int _DGEMV_(const char* transpose, ftnlen trsize,
81  const integer* m, const integer* n, const doublereal* alpha,
82  const doublereal* a, const integer* lda, const doublereal* x,
83  const integer* incX, const doublereal* beta, doublereal* y,
84  const integer* incY);
85 #endif
86 
87  int _DGETRF_(const integer* m, const integer* n,
88  doublereal* a, integer* lda, integer* ipiv,
89  integer* info);
90 
91 #ifdef LAPACK_FTN_STRING_LEN_AT_END
92  int _DGETRS_(const char* transpose, const integer* n,
93  const integer* nrhs, doublereal* a, const integer* lda,
94  integer* ipiv, doublereal* b, const integer* ldb,
95  integer* info, ftnlen trsize);
96 #else
97  int _DGETRS_(const char* transpose, ftnlen trsize, const integer* n,
98  const integer* nrhs, const doublereal* a, const integer* lda,
99  integer* ipiv, doublereal* b, const integer* ldb, integer* info);
100 #endif
101 
102  int _DGETRI_(const integer* n, doublereal* a, const integer* lda,
103  integer* ipiv, doublereal* work, integer* lwork, integer* info);
104 
105  int _DGELSS_(const integer* m, const integer* n, const integer* nrhs,
106  doublereal* a, const integer* lda, doublereal* b,
107  const integer* ldb, doublereal* s, const doublereal* rcond,
108  integer* rank, doublereal* work, integer* lwork, integer* info);
109 
110  int _DGBSV_(integer* n, integer* kl, integer* ku, integer* nrhs,
111  doublereal* a, integer* lda, integer* ipiv, doublereal* b,
112  integer* ldb, integer* info);
113 
114  int _DGBTRF_(integer* m, integer* n, integer* kl, integer* ku,
115  doublereal* a, integer* lda, integer* ipiv, integer* info);
116 
117 #ifdef LAPACK_FTN_STRING_LEN_AT_END
118  int _DGBTRS_(const char* trans, integer* n, integer* kl, integer* ku,
119  integer* nrhs, doublereal* a, integer* lda, integer* ipiv,
120  doublereal* b, integer* ldb, integer* info, ftnlen trsize);
121 #else
122  int _DGBTRS_(const char* trans, ftnlen trsize,
123  integer* n, integer* kl, integer* ku,
124  integer* nrhs, doublereal* a, integer* lda, integer* ipiv,
125  doublereal* b, integer* ldb, integer* info);
126 #endif
127 
128  int _DSCAL_(integer* n, doublereal* da, doublereal* dx, integer* incx);
129 
130  int _DGEQRF_(const integer* m, const integer* n, doublereal* a, const integer* lda,
131  doublereal* tau, doublereal* work, const integer* lwork, integer* info);
132 
133 #ifdef LAPACK_FTN_STRING_LEN_AT_END
134  int _DORMQR_(const char* side, const char* trans, const integer* m, const integer* n,
135  const integer* k, doublereal* a, const integer* lda,
136  doublereal* tau, doublereal* c, const integer* ldc,
137  doublereal* work, const integer* lwork, integer* info, ftnlen sisize, ftnlen trsize);
138 #else
139  int _DORMQR_(const char* side, ftnlen sisize, const char* trans, ftnlen trsize, const integer* m,
140  const integer* n, const integer* k, doublereal* a, const integer* lda,
141  doublereal* tau,doublereal* c, const integer* ldc,
142  doublereal* work, const integer* lwork, integer* info);
143 #endif
144 
145 #ifdef LAPACK_FTN_STRING_LEN_AT_END
146  int _DTRTRS_(const char* uplo, const char* trans, const char* diag, const integer* n,
147  const integer* nrhs, doublereal* a, const integer* lda,
148  doublereal* b, const integer* ldb, integer* info,
149  ftnlen upsize, ftnlen trsize, ftnlen disize);
150 #else
151  int _DTRTRS_(const char* uplo, ftnlen upsize, const char* trans, ftnlen trsize, const char* diag,
152  ftnlen disize, const integer* n, const integer* nrhs, doublereal* a, const integer* lda,
153  doublereal* b, const integer* ldb, integer* info);
154 #endif
155 
156 
157 #ifdef LAPACK_FTN_STRING_LEN_AT_END
158  int _DTRCON_(const char* norm, const char* uplo, const char* diag, const integer* n,
159  doublereal* a, const integer* lda, const doublereal* rcond,
160  doublereal* work, const integer* iwork, integer* info, ftnlen nosize,
161  ftnlen upsize, ftnlen disize);
162 #else
163  int _DTRCON_(const char* norm, ftnlen nosize, const char* uplo, ftnlen upsize, const char* diag,
164  ftnlen disize, const integer* n, doublereal* a, const integer* lda, const doublereal* rcond,
165  doublereal* work, const integer* iwork, integer* info);
166 #endif
167 
168 
169 #ifdef LAPACK_FTN_STRING_LEN_AT_END
170  int _DPOTRF_(const char* uplo, const integer* n, doublereal* a, const integer* lda, integer* info,
171  ftnlen upsize);
172 #else
173  int _DPOTRF_(const char* uplo, ftnlen upsize, const integer* n, doublereal* a, const integer* lda,
174  integer* info);
175 #endif
176 
177 #ifdef LAPACK_FTN_STRING_LEN_AT_END
178  int _DPOTRS_(const char* uplo, const integer* n, const integer* nrhs, doublereal* a, const integer* lda,
179  doublereal* b, const integer* ldb, integer* info, ftnlen upsize);
180 #else
181  int _DPOTRS_(const char* uplo, ftnlen upsize, const integer* n, const integer* nrhs, doublereal* a, const integer* lda,
182  doublereal* b, const integer* ldb, integer* info);
183 #endif
184 
185 #ifdef LAPACK_FTN_STRING_LEN_AT_END
186  int _DGECON_(const char* norm, const integer* n, doublereal* a, const integer* lda,
187  const doublereal* rnorm, const doublereal* rcond,
188  doublereal* work, const integer* iwork, integer* info, ftnlen nosize);
189 #else
190  int _DGECON_(const char* norm, ftnlen nosize, const integer* n, doublereal* a, const integer* lda,
191  const doublereal* rnorm, const doublereal* rcond,
192  doublereal* work, const integer* iwork, integer* info);
193 #endif
194 
195 #ifdef LAPACK_FTN_STRING_LEN_AT_END
196  int _DGBCON_(const char* norm, const integer* n, integer* kl, integer* ku, doublereal* ab, const integer* ldab,
197  const integer* ipiv, const doublereal* anorm, const doublereal* rcond,
198  doublereal* work, const integer* iwork, integer* info, ftnlen nosize);
199 #else
200  int _DGBCON_(const char* norm, ftnlen nosize, const integer* n, integer* kl, integer* ku, doublereal* ab, const integer* ldab,
201  const integer* ipiv, const doublereal* anorm, const doublereal* rcond,
202  doublereal* work, const integer* iwork, integer* info);
203 #endif
204 
205 #ifdef LAPACK_FTN_STRING_LEN_AT_END
206  doublereal _DLANGE_(const char* norm, const integer* m, const integer* n, doublereal* a, const integer* lda,
207  doublereal* work, ftnlen nosize);
208 #else
209  doublereal _DLANGE_(const char* norm, ftnlen nosize, const integer* m, const integer* n, doublereal* a, const integer* lda,
210  doublereal* work);
211 #endif
212 }
213 
214 namespace Cantera
215 {
216 
217 inline void ct_dgemv(ctlapack::storage_t storage,
218  ctlapack::transpose_t trans,
219  int m, int n, doublereal alpha, const doublereal* a, int lda,
220  const doublereal* x, int incX, doublereal beta,
221  doublereal* y, int incY)
222 {
223  integer f_m = m, f_n = n, f_lda = lda, f_incX = incX, f_incY = incY;
224  doublereal f_alpha = alpha, f_beta = beta;
225  ftnlen trsize = 1;
226 #ifdef LAPACK_FTN_STRING_LEN_AT_END
227  _DGEMV_(&no_yes[trans], &f_m, &f_n, &f_alpha, a,
228  &f_lda, x, &f_incX, &f_beta, y, &f_incY, trsize);
229 #else
230  _DGEMV_(&no_yes[trans], trsize, &f_m, &f_n, &f_alpha, a,
231  &f_lda, x, &f_incX, &f_beta, y, &f_incY);
232 #endif
233 }
234 
235 inline void ct_dgbsv(int n, int kl, int ku, int nrhs,
236  doublereal* a, int lda, integer* ipiv, doublereal* b, int ldb,
237  int& info)
238 {
239  integer f_n = n, f_kl = kl, f_ku = ku, f_nrhs = nrhs, f_lda = lda,
240  f_ldb = ldb, f_info = 0;
241  _DGBSV_(&f_n, &f_kl, &f_ku, &f_nrhs, a, &f_lda, ipiv,
242  b, &f_ldb, &f_info);
243  info = f_info;
244 }
245 
246 inline void ct_dgelss(size_t m, size_t n, size_t nrhs, doublereal* a,
247  size_t lda, doublereal* b, size_t ldb, doublereal* s,
248  doublereal rcond, size_t& rank, doublereal* work,
249  int& lwork, int& info)
250 {
251  integer f_m = static_cast<integer>(m);
252  integer f_n = static_cast<integer>(n);
253  integer f_nrhs = static_cast<integer>(nrhs);
254  integer f_lda = static_cast<integer>(lda);
255  integer f_ldb = static_cast<integer>(ldb);
256  integer f_lwork = static_cast<integer>(lwork);
257  integer f_rank, f_info;
258 
259  _DGELSS_(&f_m, &f_n, &f_nrhs, a, &f_lda, b, &f_ldb, s, &rcond, &f_rank,
260  work, &f_lwork, &f_info);
261 
262  info = static_cast<int>(f_info);
263  rank = static_cast<size_t>(f_rank);
264  lwork = static_cast<int>(f_lwork);
265 }
266 
267 inline void ct_dgbtrf(size_t m, size_t n, size_t kl, size_t ku,
268  doublereal* a, size_t lda, integer* ipiv, int& info)
269 {
270  integer f_m = (int) m;
271  integer f_n = (int) n;
272  integer f_kl = (int) kl;
273  integer f_ku = (int) ku;
274  integer f_lda = (int) lda;
275  integer f_info = 0;
276  _DGBTRF_(&f_m, &f_n, &f_kl, &f_ku, a, &f_lda, ipiv, &f_info);
277  info = f_info;
278 }
279 
280 inline void ct_dgbtrs(ctlapack::transpose_t trans, size_t n,
281  size_t kl, size_t ku, size_t nrhs, doublereal* a, size_t lda,
282  integer* ipiv, doublereal* b, size_t ldb, int& info)
283 {
284  integer f_n = (int) n;
285  integer f_kl = (int) kl;
286  integer f_ku = (int) ku;
287  integer f_nrhs = (int) nrhs;
288  integer f_lda = (int) lda;
289  integer f_ldb = (int) ldb;
290  integer f_info = 0;
291  char tr = no_yes[trans];
292  ftnlen trsize = 1;
293 #ifdef LAPACK_FTN_STRING_LEN_AT_END
294  _DGBTRS_(&tr, &f_n, &f_kl, &f_ku, &f_nrhs, a, &f_lda, ipiv,
295  b, &f_ldb, &f_info, trsize);
296 #else
297  _DGBTRS_(&tr, trsize, &f_n, &f_kl, &f_ku, &f_nrhs, a, &f_lda, ipiv,
298  b, &f_ldb, &f_info);
299 #endif
300  info = f_info;
301 }
302 
303 inline void ct_dgetrf(size_t m, size_t n,
304  doublereal* a, size_t lda, integer* ipiv, int& info)
305 {
306  integer mm = (int) m;
307  integer nn = (int) n;
308  integer ldaa = (int) lda;
309  integer infoo = 0;
310  _DGETRF_(&mm, &nn, a, &ldaa, ipiv, &infoo);
311  info = infoo;
312 }
313 
314 inline void ct_dgetrs(ctlapack::transpose_t trans, size_t n,
315  size_t nrhs, doublereal* a, size_t lda,
316  integer* ipiv, doublereal* b, size_t ldb, int& info)
317 {
318  integer f_n = (int) n;
319  integer f_lda = (int) lda;
320  integer f_nrhs = (int) nrhs;
321  integer f_ldb = (int) ldb;
322  integer f_info = 0;
323  char tr = no_yes[trans];
324  ftnlen trsize = 1;
325 #ifdef LAPACK_FTN_STRING_LEN_AT_END
326  _DGETRS_(&tr, &f_n, &f_nrhs, a, &f_lda, ipiv, b, &f_ldb, &f_info, trsize);
327 #else
328  _DGETRS_(&tr, trsize, &f_n, &f_nrhs, a, &f_lda, ipiv, b, &f_ldb, &f_info);
329 #endif
330  info = f_info;
331 }
332 
333 inline void ct_dgetri(int n, doublereal* a, int lda, integer* ipiv,
334  doublereal* work, int lwork, int& info)
335 {
336  integer f_n = n, f_lda = lda, f_lwork = lwork, f_info = 0;
337  _DGETRI_(&f_n, a, &f_lda, ipiv, work, &f_lwork, &f_info);
338 }
339 
340 inline void ct_dscal(int n, doublereal da, doublereal* dx, int incx)
341 {
342  integer f_n = n, f_incx = incx;
343  _DSCAL_(&f_n, &da, dx, &f_incx);
344 }
345 
346 inline void ct_dgeqrf(size_t m, size_t n, doublereal* a, size_t lda, doublereal* tau,
347  doublereal* work, size_t lwork, int& info)
348 {
349  integer f_m = static_cast<integer>(m);
350  integer f_n = static_cast<integer>(n);
351  integer f_lda = static_cast<integer>(lda);
352  integer f_lwork = static_cast<integer>(lwork);
353  integer f_info = 0;
354  _DGEQRF_(&f_m, &f_n, a, &f_lda, tau, work, &f_lwork, &f_info);
355  info = f_info;
356 }
357 
358 inline void ct_dormqr(ctlapack::side_t rlside, ctlapack::transpose_t trans, size_t m,
359  size_t n, size_t k, doublereal* a, size_t lda, doublereal* tau, doublereal* c, size_t ldc,
360  doublereal* work, size_t lwork, int& info)
361 {
362  char side = left_right[rlside];
363  char tr = no_yes[trans];
364  integer f_m = static_cast<integer>(m);
365  integer f_n = static_cast<integer>(n);
366  integer f_k = static_cast<integer>(k);
367  integer f_lwork = static_cast<integer>(lwork);
368  integer f_lda = static_cast<integer>(lda);
369  integer f_ldc = static_cast<integer>(ldc);
370  integer f_info = 0;
371  ftnlen trsize = 1;
372 #ifdef LAPACK_FTN_STRING_LEN_AT_END
373  _DORMQR_(&side, &tr, &f_m, &f_n, &f_k, a, &f_lda, tau, c, &f_ldc, work, &f_lwork, &f_info, trsize, trsize);
374 #else
375  _DORMQR_(&side, trsize, &tr, trsize, &f_m, &f_n, &f_k, a, &f_lda, tau, c, &f_ldc, work, &f_lwork, &f_info);
376 #endif
377  info = f_info;
378 }
379 
380 inline void ct_dtrtrs(ctlapack::upperlower_t uplot, ctlapack::transpose_t trans, const char* diag,
381  size_t n, size_t nrhs, doublereal* a, size_t lda, doublereal* b, size_t ldb, int& info)
382 {
383  char uplo = upper_lower[uplot];
384  char tr = no_yes[trans];
385  char dd = 'N';
386  if (diag) {
387  dd = diag[0];
388  }
389  integer f_n = static_cast<integer>(n);
390  integer f_nrhs = static_cast<integer>(nrhs);
391  integer f_lda = static_cast<integer>(lda);
392  integer f_ldb = static_cast<integer>(ldb);
393  integer f_info = 0;
394  ftnlen trsize = 1;
395 #ifdef LAPACK_FTN_STRING_LEN_AT_END
396  _DTRTRS_(&uplo, &tr, &dd, &f_n, &f_nrhs, a, &f_lda, b, &f_ldb, &f_info, trsize, trsize, trsize);
397 #else
398  _DTRTRS_(&uplo, trsize, &tr, trsize, &dd, trsize, &f_n, &f_nrhs, a, &f_lda, b, &f_ldb, &f_info);
399 #endif
400  info = f_info;
401 }
402 
403 /**
404  * - `work` must be dimensioned equal to greater than 3N
405  * - `iwork` must be dimensioned equal to or greater than N
406  */
407 inline doublereal ct_dtrcon(const char* norm, ctlapack::upperlower_t uplot, const char* diag,
408  size_t n, doublereal* a, size_t lda, doublereal* work, int* iwork, int& info)
409 {
410  char uplo = upper_lower[uplot];
411  char dd = 'N';
412  if (diag) {
413  dd = diag[0];
414  }
415  char nn = '1';
416  if (norm) {
417  nn = norm[0];
418  }
419  integer f_n = static_cast<integer>(n);
420  integer f_lda = static_cast<integer>(lda);
421  integer f_info = 0;
422  doublereal rcond = 0.0;
423  ftnlen trsize = 1;
424 #ifdef LAPACK_FTN_STRING_LEN_AT_END
425  _DTRCON_(&nn, &uplo, &dd, &f_n, a, &f_lda, &rcond, work, iwork, &f_info, trsize, trsize, trsize);
426 #else
427  _DTRCON_(&nn, trsize, &uplo, trsize, &dd, trsize, &f_n, a, &f_lda, &rcond, work, iwork, &f_info);
428 #endif
429  info = f_info;
430  return rcond;
431 }
432 
433 inline void ct_dpotrf(ctlapack::upperlower_t uplot, size_t n, doublereal* a, size_t lda, int& info)
434 {
435  char uplo = upper_lower[uplot];
436  integer f_n = static_cast<integer>(n);
437  integer f_lda = static_cast<integer>(lda);
438  integer f_info = 0;
439  ftnlen trsize = 1;
440 #ifdef LAPACK_FTN_STRING_LEN_AT_END
441  _DPOTRF_(&uplo, &f_n, a, &f_lda, &f_info, trsize);
442 #else
443  _DPOTRF_(&uplo, trsize, &f_n, a, &f_lda, &f_info);
444 #endif
445  info = f_info;
446  return;
447 }
448 
449 inline void ct_dpotrs(ctlapack::upperlower_t uplot, size_t n, size_t nrhs, doublereal* a, size_t lda,
450  doublereal* b, size_t ldb, int& info)
451 {
452  char uplo = upper_lower[uplot];
453  integer f_n = static_cast<integer>(n);
454  integer f_nrhs = static_cast<integer>(nrhs);
455  integer f_lda = static_cast<integer>(lda);
456  integer f_ldb = static_cast<integer>(ldb);
457  integer f_info = 0;
458  ftnlen trsize = 1;
459 #ifdef LAPACK_FTN_STRING_LEN_AT_END
460  _DPOTRS_(&uplo, &f_n, &f_nrhs, a, &f_lda, b, &f_ldb, &f_info, trsize);
461 #else
462  _DPOTRS_(&uplo, trsize, &f_n, &f_nrhs, a, &f_lda, b, &f_ldb, &f_info);
463 #endif
464  info = f_info;
465  return;
466 }
467 
468 inline doublereal ct_dgecon(const char norm, size_t n, doublereal* a, size_t lda, doublereal anorm,
469  doublereal* work, int* iwork, int& info)
470 {
471  char cnorm = '1';
472  if (norm) {
473  cnorm = norm;
474  }
475  integer f_n = static_cast<integer>(n);
476  integer f_lda = static_cast<integer>(lda);
477  integer f_info = 0;
478  doublereal rcond = 0.0;
479  ftnlen trsize = 1;
480 #ifdef LAPACK_FTN_STRING_LEN_AT_END
481  _DGECON_(&cnorm, &f_n, a, &f_lda, &anorm, &rcond, work, iwork, &f_info, trsize);
482 #else
483  _DGECON_(&cnorm, trsize, &f_n, a, &f_lda, &anorm, &rcond, work, iwork, &f_info);
484 #endif
485  info = f_info;
486  return rcond;
487 }
488 
489 inline doublereal ct_dgbcon(const char norm, size_t n, size_t kl, size_t ku,
490  doublereal* a, size_t ldab, int* ipiv, doublereal anorm,
491  doublereal* work, int* iwork, int& info)
492 {
493  char cnorm = '1';
494  if (norm) {
495  cnorm = norm;
496  }
497  integer f_n = static_cast<integer>(n);
498  integer f_kl = static_cast<integer>(kl);
499  integer f_ku = static_cast<integer>(ku);
500  integer f_ldab = static_cast<integer>(ldab);
501  integer f_info = 0;
502  doublereal rcond = 0.0;
503  ftnlen trsize = 1;
504 #ifdef LAPACK_FTN_STRING_LEN_AT_END
505  _DGBCON_(&cnorm, &f_n, &f_kl, &f_ku, a, &f_ldab, ipiv, &anorm, &rcond, work, iwork, &f_info, trsize);
506 #else
507  _DGBCON_(&cnorm, trsize, &f_n, &f_kl, &f_ku, a, &f_ldab, ipiv, &anorm, &rcond, work, iwork, &f_info);
508 #endif
509  info = f_info;
510  return rcond;
511 }
512 
513 inline doublereal ct_dlange(const char norm, size_t m, size_t n, doublereal* a, size_t lda,
514  doublereal* work)
515 {
516  char cnorm = '1';
517  if (norm) {
518  cnorm = norm;
519  }
520  integer f_m = static_cast<integer>(m);
521  integer f_n = static_cast<integer>(n);
522  integer f_lda = static_cast<integer>(lda);
523  doublereal anorm;
524  ftnlen trsize = 1;
525 #ifdef LAPACK_FTN_STRING_LEN_AT_END
526  anorm = _DLANGE_(&cnorm, &f_m, &f_n, a, &f_lda, work, trsize);
527 #else
528  anorm = _DLANGE_(&cnorm, trsize, &f_m, &f_n, a, &f_lda, work);
529 #endif
530  return anorm;
531 }
532 
533 }
534 
535 #endif
This file contains definitions of constants, types and terms that are used in internal routines and a...
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:564
doublereal ct_dtrcon(const char *norm, ctlapack::upperlower_t uplot, const char *diag, size_t n, doublereal *a, size_t lda, doublereal *work, int *iwork, int &info)
Definition: ctlapack.h:407