QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_c_interface.F90
Go to the documentation of this file.
1 !! ##############################################################################################
2 !!
3 !! Copyright 2012 CNRS, INPT
4 !!
5 !! This file is part of qr_mumps.
6 !!
7 !! qr_mumps is free software: you can redistribute it and/or modify
8 !! it under the terms of the GNU Lesser General Public License as
9 !! published by the Free Software Foundation, either version 3 of
10 !! the License, or (at your option) any later version.
11 !!
12 !! qr_mumps is distributed in the hope that it will be useful,
13 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
14 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 !! GNU Lesser General Public License for more details.
16 !!
17 !! You can find a copy of the GNU Lesser General Public License
18 !! in the qr_mumps/doc directory.
19 !!
20 !! ##############################################################################################
21 
22 
23 !! ##############################################################################################
33 
34 
35 #include "qrm_common.h"
36 
39  use iso_c_binding
40  use _qrm_spmat_mod
41 
45  type, bind(c) :: _qrm_spmat_type_c
48  type(c_ptr) :: irn
49 
51  type(c_ptr) :: jcn
52 
54  type(c_ptr) :: val
55 
56  integer(c_int) :: m
57 
58  integer(c_int) :: n
59 
60  integer(c_int) :: nz
61 
63  type(c_ptr) :: cperm_in, rperm, cperm
64 
65  integer(c_int) :: icntl(20)
66 
67  real(c_double) :: rcntl(10)
68 
69  integer(c_long) :: gstats(10)
70 
72  integer(c_int) :: h
73  type(c_ptr) :: mat_ptr
74  end type _qrm_spmat_type_c
75 
76  type inst
77  type(_qrm_spmat_type), pointer :: fmat => null()
78  end type inst
79 
80  integer, parameter :: max_inst=10
81  type(inst), save :: spmat_array(max_inst)
82  integer, save :: ninst=0
83 
84 contains
85 
88  subroutine _qrm_spmat_init_c(qrm_spmat_c) bind(c)
89  use _qrm_spmat_mod
90  implicit none
91  type(_qrm_spmat_type_c) :: qrm_spmat_c
92  type(_qrm_spmat_type), pointer :: fmat
93  integer :: i
94 
95  allocate(fmat)
96  call _qrm_spmat_init(fmat)
97 
98  qrm_spmat_c%icntl = fmat%icntl
99  qrm_spmat_c%rcntl = fmat%rcntl
100  qrm_spmat_c%gstats = fmat%gstats
101 
102  qrm_spmat_c%mat_ptr = c_loc(fmat)
103 
104  nullify(fmat)
105 
106  return
107 
108  end subroutine _qrm_spmat_init_c
109 
110 
113  subroutine _qrm_spmat_destroy_c(qrm_spmat_c) bind(c)
114  use _qrm_spmat_mod
115  implicit none
116  type(_qrm_spmat_type_c) :: qrm_spmat_c
117 
118  type(_qrm_spmat_type), pointer :: fmat
119  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
120 
121  call _qrm_spmat_destroy(fmat, all=.false.)
122 
123  qrm_spmat_c%icntl = fmat%icntl
124  qrm_spmat_c%rcntl = fmat%rcntl
125  qrm_spmat_c%gstats = fmat%gstats
126 
127  qrm_spmat_c%mat_ptr = c_null_ptr
128  return
129 
130  end subroutine _qrm_spmat_destroy_c
131 
132 
135  subroutine _qrm_analyse_c(qrm_spmat_c, transp) bind(c)
137  implicit none
138  type(_qrm_spmat_type_c) :: qrm_spmat_c
139  character(kind=c_char), value :: transp
140 
141  type(_qrm_spmat_type), pointer :: fmat
142  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
143 
144  fmat%m = qrm_spmat_c%m
145  fmat%n = qrm_spmat_c%n
146  fmat%nz = qrm_spmat_c%nz
147 
148  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
149  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
150  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
151  call c_f_pointer(qrm_spmat_c%cperm_in, fmat%cperm_in,(/qrm_spmat_c%n/))
152 
153  fmat%icntl = qrm_spmat_c%icntl
154  fmat%rcntl = qrm_spmat_c%rcntl
155 
156  call _qrm_analyse(fmat, transp)
157 
158  qrm_spmat_c%gstats = fmat%gstats
159  qrm_spmat_c%cperm = c_loc(fmat%adata%cperm(1))
160  qrm_spmat_c%rperm = c_loc(fmat%adata%rperm(1))
161 
162  return
163  end subroutine _qrm_analyse_c
164 
167  subroutine _qrm_factorize_c(qrm_spmat_c, transp) bind(c)
169  implicit none
170  type(_qrm_spmat_type_c) :: qrm_spmat_c
171  character(kind=c_char), value :: transp
172 
173  type(_qrm_spmat_type), pointer :: fmat
174  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
175 
176  fmat%m = qrm_spmat_c%m
177  fmat%n = qrm_spmat_c%n
178  fmat%nz = qrm_spmat_c%nz
179 
180  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
181  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
182  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
183 
184  fmat%icntl = qrm_spmat_c%icntl
185  fmat%rcntl = qrm_spmat_c%rcntl
186 
187  call _qrm_factorize(fmat, transp)
188 
189  qrm_spmat_c%gstats = fmat%gstats
190 
191  return
192  end subroutine _qrm_factorize_c
193 
194 
197  subroutine _qrm_get_r_c(qrm_spmat_c, r_c) bind(c)
199  implicit none
200  type(_qrm_spmat_type_c) :: qrm_spmat_c, r_c
201 
202  type(_qrm_spmat_type), pointer :: fmat, rmat
203 
204  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
205  call c_f_pointer(r_c%mat_ptr, rmat)
206 
207  call _qrm_get_r(fmat, rmat)
208 
209  r_c%m = rmat%m
210  r_c%n = rmat%n
211  r_c%nz = rmat%nz
212 
213  r_c%cperm = c_loc(rmat%adata%cperm(1))
214  r_c%rperm = c_loc(rmat%adata%rperm(1))
215 
216  r_c%irn = c_loc(rmat%irn(1))
217  r_c%jcn = c_loc(rmat%jcn(1))
218  r_c%val = c_loc(rmat%val(1))
219 
220  return
221  end subroutine _qrm_get_r_c
222 
223 
226  subroutine _qrm_solve_c(qrm_spmat_c, transp, b, x, nrhs) bind(c)
227  use _qrm_solve_mod
228  implicit none
229  type(_qrm_spmat_type_c) :: qrm_spmat_c
230  character(kind=c_char), value :: transp
231  type(c_ptr), value :: b, x
232  integer(c_int), value :: nrhs
233 
234  _qrm_data, pointer :: ib(:,:), ix(:,:)
235 
236  type(_qrm_spmat_type), pointer :: fmat
237  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
238 
239  fmat%m = qrm_spmat_c%m
240  fmat%n = qrm_spmat_c%n
241  fmat%nz = qrm_spmat_c%nz
242 
243  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
244  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
245  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
246 
247  call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
248  call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
249 
250  fmat%icntl = qrm_spmat_c%icntl
251  fmat%rcntl = qrm_spmat_c%rcntl
252 
253  call _qrm_solve(fmat, transp, ib, ix)
254 
255  qrm_spmat_c%gstats = fmat%gstats
256 
257  return
258  end subroutine _qrm_solve_c
259 
260 
263  subroutine _qrm_apply_c(qrm_spmat_c, transp, b, nrhs) bind(c)
264  use _qrm_solve_mod
265  implicit none
266  type(_qrm_spmat_type_c) :: qrm_spmat_c
267  character(kind=c_char), value :: transp
268  type(c_ptr), value :: b
269  integer(c_int), value :: nrhs
270 
271  _qrm_data, pointer :: ib(:,:)
272 
273  type(_qrm_spmat_type), pointer :: fmat
274  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
275 
276  fmat%m = qrm_spmat_c%m
277  fmat%n = qrm_spmat_c%n
278  fmat%nz = qrm_spmat_c%nz
279 
280  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
281  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
282  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
283 
284  if(transp.eq.'t') then
285  call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
286  else if(transp.eq.'n') then
287  call c_f_pointer(b, ib,(/qrm_spmat_c%n, nrhs/))
288  end if
289 
290  fmat%icntl = qrm_spmat_c%icntl
291  fmat%rcntl = qrm_spmat_c%rcntl
292 
293  call _qrm_apply(fmat, transp, ib)
294 
295  qrm_spmat_c%gstats = fmat%gstats
296 
297  return
298  end subroutine _qrm_apply_c
299 
300 
301 
304  subroutine _qrm_matmul_c(qrm_spmat_c, transp, alpha, x, beta, y, nrhs) bind(c)
305  use _qrm_utils_mod
306  implicit none
307  type(_qrm_spmat_type_c) :: qrm_spmat_c
308  character(kind=c_char), value :: transp
309  _qrm_data_fc, value :: alpha, beta
310  type(c_ptr), value :: x, y
311  integer(c_int), value :: nrhs
312 
313  _qrm_data, pointer :: ix(:,:), iy(:,:)
314 
315  type(_qrm_spmat_type), pointer :: fmat
316  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
317 
318  fmat%m = qrm_spmat_c%m
319  fmat%n = qrm_spmat_c%n
320  fmat%nz = qrm_spmat_c%nz
321 
322  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
323  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
324  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
325 
326  call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
327  call c_f_pointer(y, iy,(/qrm_spmat_c%m, nrhs/))
328 
329  fmat%icntl = qrm_spmat_c%icntl
330  fmat%rcntl = qrm_spmat_c%rcntl
331 
332  call _qrm_matmul(fmat, transp, alpha, ix, beta, iy)
333 
334  qrm_spmat_c%gstats = fmat%gstats
335 
336  return
337  end subroutine _qrm_matmul_c
338 
339 
342  subroutine _qrm_matnrm_c(qrm_spmat_c, ntype, nrm) bind(c)
343  use _qrm_utils_mod
344  implicit none
345  type(_qrm_spmat_type_c) :: qrm_spmat_c
346  character(kind=c_char), value :: ntype
347  _qrm_real_fc :: nrm
348 
349  type(_qrm_spmat_type), pointer :: fmat
350  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
351 
352  fmat%m = qrm_spmat_c%m
353  fmat%n = qrm_spmat_c%n
354  fmat%nz = qrm_spmat_c%nz
355 
356  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
357  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
358  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
359 
360  fmat%icntl = qrm_spmat_c%icntl
361  fmat%rcntl = qrm_spmat_c%rcntl
362 
363  call _qrm_matnrm(fmat, ntype, nrm)
364 
365  return
366  end subroutine _qrm_matnrm_c
367 
370  subroutine _qrm_vecnrm_c(x, n, nrhs, ntype, nrm) bind(c)
371  use _qrm_utils_mod
372  implicit none
373  type(c_ptr), value :: x
374  integer(c_int), value :: nrhs
375  integer(c_int), value :: n
376  character(kind=c_char), value :: ntype
377  type(c_ptr), value :: nrm
378  _qrm_data, pointer :: ix(:,:)
379  _qrm_real, pointer :: inrm(:)
380 
381  call c_f_pointer(x, ix,(/n,nrhs/))
382  call c_f_pointer(nrm, inrm,(/nrhs/))
383 
384  call _qrm_vecnrm(ix, n, ntype, inrm)
385 
386  return
387  end subroutine _qrm_vecnrm_c
388 
389 
392  subroutine _qrm_least_squares_c(qrm_spmat_c, b, x, nrhs) bind(c)
393  use _qrm_spmat_mod
394  use _qrm_methods_mod
395  implicit none
396  type(_qrm_spmat_type_c) :: qrm_spmat_c
397  type(c_ptr), value :: b, x
398  integer(c_int), value :: nrhs
399 
400  _qrm_data, pointer :: ib(:,:), ix(:,:)
401 
402  type(_qrm_spmat_type), pointer :: fmat
403  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
404 
405  fmat%m = qrm_spmat_c%m
406  fmat%n = qrm_spmat_c%n
407  fmat%nz = qrm_spmat_c%nz
408 
409  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
410  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
411  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
412 
413  call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
414  call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
415 
416  fmat%icntl = qrm_spmat_c%icntl
417  fmat%rcntl = qrm_spmat_c%rcntl
418 
419  call qrm_least_squares(fmat, ib, ix)
420 
421  qrm_spmat_c%gstats = fmat%gstats
422 
423  return
424  end subroutine _qrm_least_squares_c
425 
428  subroutine _qrm_min_norm_c(qrm_spmat_c, b, x, nrhs) bind(c)
429  use _qrm_spmat_mod
430  use _qrm_methods_mod
431  implicit none
432  type(_qrm_spmat_type_c) :: qrm_spmat_c
433  type(c_ptr), value :: b, x
434  integer(c_int), value :: nrhs
435 
436  _qrm_data, pointer :: ib(:,:), ix(:,:)
437 
438  type(_qrm_spmat_type), pointer :: fmat
439  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
440 
441  fmat%m = qrm_spmat_c%m
442  fmat%n = qrm_spmat_c%n
443  fmat%nz = qrm_spmat_c%nz
444 
445  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
446  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
447  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
448 
449  call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
450  call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
451 
452  fmat%icntl = qrm_spmat_c%icntl
453  fmat%rcntl = qrm_spmat_c%rcntl
454 
455  call qrm_min_norm(fmat, ib, ix)
456 
457  qrm_spmat_c%gstats = fmat%gstats
458 
459  return
460  end subroutine _qrm_min_norm_c
461 
462 
465  subroutine _qrm_residual_norm_c(qrm_spmat_c, b, x, nrhs, nrm) bind(c)
466  use _qrm_spmat_mod
467  use _qrm_methods_mod
468  implicit none
469  type(_qrm_spmat_type_c) :: qrm_spmat_c
470  integer(c_int), value :: nrhs
471  type(c_ptr), value :: b, x
472  type(c_ptr), value :: nrm
473 
474  _qrm_data, pointer :: ib(:,:), ix(:,:)
475  _qrm_real, pointer :: inrm(:)
476 
477  type(_qrm_spmat_type), pointer :: fmat
478  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
479 
480  fmat%m = qrm_spmat_c%m
481  fmat%n = qrm_spmat_c%n
482  fmat%nz = qrm_spmat_c%nz
483 
484  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
485  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
486  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
487 
488  call c_f_pointer(b, ib,(/qrm_spmat_c%m,nrhs/))
489  call c_f_pointer(x, ix,(/qrm_spmat_c%n,nrhs/))
490  call c_f_pointer(nrm, inrm,(/nrhs/))
491 
492  fmat%icntl = qrm_spmat_c%icntl
493  fmat%rcntl = qrm_spmat_c%rcntl
494 
495  call qrm_residual_norm(fmat, ib, ix, inrm)
496 
497  qrm_spmat_c%gstats = fmat%gstats
498 
499  return
500  end subroutine _qrm_residual_norm_c
501 
504  subroutine _qrm_residual_orth_c(qrm_spmat_c, r, nrhs, nrm) bind(c)
505  use _qrm_spmat_mod
506  use _qrm_methods_mod
507  implicit none
508  type(_qrm_spmat_type_c) :: qrm_spmat_c
509  type(c_ptr), value :: r
510  integer(c_int), value :: nrhs
511  type(c_ptr), value :: nrm
512 
513  _qrm_data, pointer :: ir(:,:)
514  _qrm_real, pointer :: inrm(:)
515 
516  type(_qrm_spmat_type), pointer :: fmat
517  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
518 
519  fmat%m = qrm_spmat_c%m
520  fmat%n = qrm_spmat_c%n
521  fmat%nz = qrm_spmat_c%nz
522 
523  call c_f_pointer(qrm_spmat_c%irn, fmat%irn,(/qrm_spmat_c%nz/))
524  call c_f_pointer(qrm_spmat_c%jcn, fmat%jcn,(/qrm_spmat_c%nz/))
525  call c_f_pointer(qrm_spmat_c%val, fmat%val,(/qrm_spmat_c%nz/))
526 
527  call c_f_pointer(r, ir,(/qrm_spmat_c%m,nrhs/))
528  call c_f_pointer(nrm, inrm,(/nrhs/))
529 
530  fmat%icntl = qrm_spmat_c%icntl
531  fmat%rcntl = qrm_spmat_c%rcntl
532 
533  call qrm_residual_orth(fmat, ir, inrm)
534 
535  qrm_spmat_c%gstats = fmat%gstats
536 
537  return
538  end subroutine _qrm_residual_orth_c
539 
540 
541 
542 
545  subroutine _qrm_pseti_c(qrm_spmat_c, string, val) bind(c)
546  use _qrm_spmat_mod
547  implicit none
548  type(_qrm_spmat_type_c) :: qrm_spmat_c
549  character(kind=c_char) :: string(40)
550  integer(c_int), value :: val
551 
552  character(len=40) :: a
553 
554  type(_qrm_spmat_type), pointer :: fmat
555  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
556 
557  write(a,'(40a)')string
558 
559  call _qrm_pseti(fmat, a, val)
560 
561  qrm_spmat_c%icntl = fmat%icntl
562  qrm_spmat_c%rcntl = fmat%rcntl
563 
564  return
565 
566  end subroutine _qrm_pseti_c
567 
568 
571  subroutine _qrm_pgeti_c(qrm_spmat_c, string, val) bind(c)
572  use _qrm_spmat_mod
573  implicit none
574  type(_qrm_spmat_type_c) :: qrm_spmat_c
575  character(kind=c_char) :: string(40)
576  integer(c_int) :: val
577 
578  character(len=40) :: a
579 
580  type(_qrm_spmat_type), pointer :: fmat
581  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
582 
583  write(a,'(40a)')string
584 
585  fmat%icntl = qrm_spmat_c%icntl
586  fmat%rcntl = qrm_spmat_c%rcntl
587 
588  call _qrm_pgeti(fmat, a, val)
589 
590  return
591 
592  end subroutine _qrm_pgeti_c
593 
596  subroutine _qrm_pgetii_c(qrm_spmat_c, string, val) bind(c)
597  use _qrm_spmat_mod
598  implicit none
599  type(_qrm_spmat_type_c) :: qrm_spmat_c
600  character(kind=c_char) :: string(40)
601  integer(c_long_long) :: val
602 
603  character(len=40) :: a
604 
605  type(_qrm_spmat_type), pointer :: fmat
606  call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
607 
608  write(a,'(40a)')string
609 
610  fmat%icntl = qrm_spmat_c%icntl
611  fmat%rcntl = qrm_spmat_c%rcntl
612 
613  call _qrm_pgetii(fmat, a, val)
614 
615  return
616 
617  end subroutine _qrm_pgetii_c
618 
619 
620 end module _qrm_c_interface
621 
void _qrm_residual_orth_c(struct _qrm_spmat_type_c *qrm_spmat_c, _qrm_data_c *r, const int nrhs, _qrm_real_c *nrm)
void _qrm_solve_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char transp, _qrm_data_c *b, _qrm_data_c *x, const int nrhs)
Generic interface for the ::_qrm_solve and ::_qrm_solve1d routines.
Generic interface for the ::_qrm_vecnrm2d and ::_qrm_vecnrm1d routines.
void _qrm_vecnrm_c(const _qrm_data_c *x, const int n, const int nrhs, const char ntype, _qrm_real_c *nrm)
void _qrm_pgeti_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char *string, int *val)
void _qrm_get_r_c(struct _qrm_spmat_type_c *qrm_spmat_c, struct _qrm_spmat_type_c *r)
This module contains all the generic interfaces for the typed routines in the factorization phase...
subroutine _qrm_get_r(qrm_mat, r)
subroutine _qrm_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
This module contains the generic interfaces for all the analysis routines.
void _qrm_analyse_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char transp)
subroutine _qrm_pseti(qrm_spmat, string, ival)
This subroutine is meant to set the integer control parameters.
void _qrm_spmat_init_c(struct _qrm_spmat_type_c *qrm_spmat_c)
This module contains generic interfaces for a number of auxiliary tools.
void _qrm_min_norm_c(struct _qrm_spmat_type_c *qrm_spmat_c, _qrm_data_c *b, _qrm_data_c *x, const int nrhs)
subroutine _qrm_matnrm(qrm_mat, ntype, nrm)
This subroutine computes the matrix norm. The return value is a real scalar.
Definition: qrm_matnrm.F90:47
subroutine _qrm_analyse(qrm_mat, transp)
This is the driver routine for the analysis phase.
Definition: qrm_analyse.F90:64
subroutine _qrm_pgeti(qrm_spmat, string, ival)
Gets the values of an integer control parameter. This is the dual of the ::_qrm_pseti routine; the pa...
void _qrm_spmat_destroy_c(struct _qrm_spmat_type_c *qrm_spmat_c)
This module contains the definition of the qr_mumps C interface.
This module contains generic methods.
subroutine _qrm_spmat_init(qrm_spmat)
This subroutine initializes a qrm_spmat_type instance setting default values into the control paramet...
This module contains all the interfaces for the typed routines in the solve phase.
void _qrm_factorize_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char transp)
void _qrm_residual_norm_c(struct _qrm_spmat_type_c *qrm_spmat_c, _qrm_data_c *b, _qrm_data_c *x, const int nrhs, _qrm_real_c *nrm)
This type defines the data structure used to store a matrix.
void _qrm_pgetii_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char *string, long long *val)
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the ::_qrm_matmul2d and ::_qrm_matmul1d routines.
subroutine _qrm_pgetii(qrm_spmat, string, ival)
Gets the values of an integer control parameter. This is the dual of the ::_qrm_pseti routine; the pa...
Generic interface for the ::_qrm_apply and ::_qrm_apply1d routines.
void _qrm_least_squares_c(struct _qrm_spmat_type_c *qrm_spmat_c, _qrm_data_c *b, _qrm_data_c *x, const int nrhs)
subroutine _qrm_factorize(qrm_mat, transp)
This routine is the main factorization driver.
int i
Definition: secs.c:40
void _qrm_apply_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char transp, _qrm_data_c *b, const int nrhs)
void _qrm_pseti_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char *string, int val)
void _qrm_matnrm_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char ntype, _qrm_real_c *nrm)
void _qrm_matmul_c(struct _qrm_spmat_type_c *qrm_spmat_c, const char transp, const _qrm_data_c alpha, _qrm_data_c *x, const _qrm_data_c beta, _qrm_data_c *y, const int nrhs)