35 #include "qrm_common.h"
63 type(c_ptr
) :: cperm_in, rperm, cperm
65 integer(c_int) :: icntl(20)
67 real(c_double) :: rcntl(10)
69 integer(c_long) :: gstats(10)
73 type(c_ptr
) :: mat_ptr
77 type(dqrm_spmat_type
),
pointer :: fmat => null()
80 integer,
parameter :: max_inst=10
81 type(inst),
save :: spmat_array(max_inst)
82 integer,
save :: ninst=0
98 qrm_spmat_c%icntl = fmat%icntl
99 qrm_spmat_c%rcntl = fmat%rcntl
100 qrm_spmat_c%gstats = fmat%gstats
102 qrm_spmat_c%mat_ptr = c_loc(fmat)
119 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
123 qrm_spmat_c%icntl = fmat%icntl
124 qrm_spmat_c%rcntl = fmat%rcntl
125 qrm_spmat_c%gstats = fmat%gstats
127 qrm_spmat_c%mat_ptr = c_null_ptr
139 character(kind=c_char), value :: transp
141 type(dqrm_spmat_type
),
pointer :: fmat
142 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
144 fmat%m = qrm_spmat_c%m
145 fmat%n = qrm_spmat_c%n
146 fmat%nz = qrm_spmat_c%nz
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/))
153 fmat%icntl = qrm_spmat_c%icntl
154 fmat%rcntl = qrm_spmat_c%rcntl
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))
171 character(kind=c_char), value :: transp
173 type(dqrm_spmat_type
),
pointer :: fmat
174 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
176 fmat%m = qrm_spmat_c%m
177 fmat%n = qrm_spmat_c%n
178 fmat%nz = qrm_spmat_c%nz
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/))
184 fmat%icntl = qrm_spmat_c%icntl
185 fmat%rcntl = qrm_spmat_c%rcntl
189 qrm_spmat_c%gstats = fmat%gstats
202 type(dqrm_spmat_type
),
pointer :: fmat, rmat
204 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
205 call c_f_pointer(r_c%mat_ptr, rmat)
213 r_c%cperm = c_loc(rmat%adata%cperm(1))
214 r_c%rperm = c_loc(rmat%adata%rperm(1))
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))
230 character(kind=c_char), value :: transp
231 type(c_ptr
), value :: b, x
232 integer(c_int), value :: nrhs
234 real(kind(1.d0)),
pointer :: ib(:,:), ix(:,:)
236 type(dqrm_spmat_type
),
pointer :: fmat
237 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
239 fmat%m = qrm_spmat_c%m
240 fmat%n = qrm_spmat_c%n
241 fmat%nz = qrm_spmat_c%nz
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/))
247 call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
248 call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
250 fmat%icntl = qrm_spmat_c%icntl
251 fmat%rcntl = qrm_spmat_c%rcntl
255 qrm_spmat_c%gstats = fmat%gstats
267 character(kind=c_char), value :: transp
268 type(c_ptr
), value :: b
269 integer(c_int), value :: nrhs
271 real(kind(1.d0)),
pointer :: ib(:,:)
273 type(dqrm_spmat_type
),
pointer :: fmat
274 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
276 fmat%m = qrm_spmat_c%m
277 fmat%n = qrm_spmat_c%n
278 fmat%nz = qrm_spmat_c%nz
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/))
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/))
290 fmat%icntl = qrm_spmat_c%icntl
291 fmat%rcntl = qrm_spmat_c%rcntl
295 qrm_spmat_c%gstats = fmat%gstats
304 subroutine dqrm_matmul_c(qrm_spmat_c, transp, alpha, x, beta, y, nrhs) bind(c)
308 character(kind=c_char), value :: transp
309 real(c_double), value :: alpha, beta
310 type(c_ptr
), value :: x, y
311 integer(c_int), value :: nrhs
313 real(kind(1.d0)),
pointer :: ix(:,:), iy(:,:)
315 type(dqrm_spmat_type
),
pointer :: fmat
316 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
318 fmat%m = qrm_spmat_c%m
319 fmat%n = qrm_spmat_c%n
320 fmat%nz = qrm_spmat_c%nz
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/))
326 call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
327 call c_f_pointer(y, iy,(/qrm_spmat_c%m, nrhs/))
329 fmat%icntl = qrm_spmat_c%icntl
330 fmat%rcntl = qrm_spmat_c%rcntl
332 call
dqrm_matmul(fmat, transp, alpha, ix, beta, iy)
334 qrm_spmat_c%gstats = fmat%gstats
346 character(kind=c_char), value :: ntype
347 real(c_double) :: nrm
349 type(dqrm_spmat_type
),
pointer :: fmat
350 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
352 fmat%m = qrm_spmat_c%m
353 fmat%n = qrm_spmat_c%n
354 fmat%nz = qrm_spmat_c%nz
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/))
360 fmat%icntl = qrm_spmat_c%icntl
361 fmat%rcntl = qrm_spmat_c%rcntl
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 real(kind(1.d0)),
pointer :: ix(:,:)
379 real(kind(1.d0)),
pointer :: inrm(:)
381 call c_f_pointer(x, ix,(/n,nrhs/))
382 call c_f_pointer(nrm, inrm,(/nrhs/))
397 type(c_ptr
), value :: b, x
398 integer(c_int), value :: nrhs
400 real(kind(1.d0)),
pointer :: ib(:,:), ix(:,:)
403 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
405 fmat%m = qrm_spmat_c%m
406 fmat%n = qrm_spmat_c%n
407 fmat%nz = qrm_spmat_c%nz
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/))
413 call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
414 call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
416 fmat%icntl = qrm_spmat_c%icntl
417 fmat%rcntl = qrm_spmat_c%rcntl
421 qrm_spmat_c%gstats = fmat%gstats
433 type(c_ptr
), value :: b, x
434 integer(c_int), value :: nrhs
436 real(kind(1.d0)),
pointer :: ib(:,:), ix(:,:)
439 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
441 fmat%m = qrm_spmat_c%m
442 fmat%n = qrm_spmat_c%n
443 fmat%nz = qrm_spmat_c%nz
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/))
449 call c_f_pointer(b, ib,(/qrm_spmat_c%m, nrhs/))
450 call c_f_pointer(x, ix,(/qrm_spmat_c%n, nrhs/))
452 fmat%icntl = qrm_spmat_c%icntl
453 fmat%rcntl = qrm_spmat_c%rcntl
457 qrm_spmat_c%gstats = fmat%gstats
470 integer(c_int), value :: nrhs
471 type(c_ptr
), value :: b, x
472 type(c_ptr
), value :: nrm
474 real(kind(1.d0)),
pointer :: ib(:,:), ix(:,:)
475 real(kind(1.d0)),
pointer :: inrm(:)
478 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
480 fmat%m = qrm_spmat_c%m
481 fmat%n = qrm_spmat_c%n
482 fmat%nz = qrm_spmat_c%nz
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/))
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/))
492 fmat%icntl = qrm_spmat_c%icntl
493 fmat%rcntl = qrm_spmat_c%rcntl
497 qrm_spmat_c%gstats = fmat%gstats
509 type(c_ptr
), value :: r
510 integer(c_int), value :: nrhs
511 type(c_ptr
), value :: nrm
513 real(kind(1.d0)),
pointer :: ir(:,:)
514 real(kind(1.d0)),
pointer :: inrm(:)
517 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
519 fmat%m = qrm_spmat_c%m
520 fmat%n = qrm_spmat_c%n
521 fmat%nz = qrm_spmat_c%nz
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/))
527 call c_f_pointer(r, ir,(/qrm_spmat_c%m,nrhs/))
528 call c_f_pointer(nrm, inrm,(/nrhs/))
530 fmat%icntl = qrm_spmat_c%icntl
531 fmat%rcntl = qrm_spmat_c%rcntl
535 qrm_spmat_c%gstats = fmat%gstats
549 character(kind=c_char) :: string(40)
550 integer(c_int), value :: val
552 character(len=40) :: a
555 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
557 write(a,
'(40a)')string
561 qrm_spmat_c%icntl = fmat%icntl
562 qrm_spmat_c%rcntl = fmat%rcntl
575 character(kind=c_char) :: string(40)
576 integer(c_int) :: val
578 character(len=40) :: a
581 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
583 write(a,
'(40a)')string
585 fmat%icntl = qrm_spmat_c%icntl
586 fmat%rcntl = qrm_spmat_c%rcntl
600 character(kind=c_char) :: string(40)
601 integer(c_long_long) :: val
603 character(len=40) :: a
606 call c_f_pointer(qrm_spmat_c%mat_ptr, fmat)
608 write(a,
'(40a)')string
610 fmat%icntl = qrm_spmat_c%icntl
611 fmat%rcntl = qrm_spmat_c%rcntl
subroutine dqrm_get_r(qrm_mat, r)
subroutine dqrm_pgetii(qrm_spmat, string, ival)
Gets the values of an integer control parameter. This is the dual of the ::dqrm_pseti routine; the pa...
This module contains all the interfaces for the typed routines in the solve phase.
void dqrm_least_squares_c(struct dqrm_spmat_type_c *qrm_spmat_c, double *b, double *x, const int nrhs)
This module contains generic methods.
void dqrm_pgeti_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char *string, int *val)
void dqrm_residual_norm_c(struct dqrm_spmat_type_c *qrm_spmat_c, double *b, double *x, const int nrhs, double *nrm)
void dqrm_pgetii_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char *string, long long *val)
Generic interface for the ::dqrm_vecnrm2d and ::dqrm_vecnrm1d routines.
void dqrm_get_r_c(struct dqrm_spmat_type_c *qrm_spmat_c, struct dqrm_spmat_type_c *r)
This module contains the definition of the qr_mumps C interface.
void dqrm_matnrm_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char ntype, double *nrm)
This module contains the definition of the basic sparse matrix type and of the associated methods...
void dqrm_solve_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char transp, double *b, double *x, const int nrhs)
This module contains generic interfaces for a number of auxiliary tools.
subroutine dqrm_pgeti(qrm_spmat, string, ival)
Gets the values of an integer control parameter. This is the dual of the ::dqrm_pseti routine; the pa...
void dqrm_vecnrm_c(const double *x, const int n, const int nrhs, const char ntype, double *nrm)
subroutine dqrm_analyse(qrm_mat, transp)
This is the driver routine for the analysis phase.
Generic interface for the ::dqrm_solve and ::dqrm_solve1d routines.
This module contains all the generic interfaces for the typed routines in the factorization phase...
void dqrm_min_norm_c(struct dqrm_spmat_type_c *qrm_spmat_c, double *b, double *x, const int nrhs)
Generic interface for the ::dqrm_matmul2d and ::dqrm_matmul1d routines.
subroutine dqrm_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
subroutine dqrm_pseti(qrm_spmat, string, ival)
This subroutine is meant to set the integer control parameters.
void dqrm_spmat_init_c(struct dqrm_spmat_type_c *qrm_spmat_c)
void dqrm_factorize_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char transp)
void dqrm_pseti_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char *string, int val)
subroutine dqrm_spmat_init(qrm_spmat)
This subroutine initializes a qrm_spmat_type instance setting default values into the control paramet...
Generic interface for the ::dqrm_apply and ::dqrm_apply1d routines.
void dqrm_apply_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char transp, double *b, const int nrhs)
This type defines the data structure used to store a matrix.
void dqrm_analyse_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char transp)
void dqrm_matmul_c(struct dqrm_spmat_type_c *qrm_spmat_c, const char transp, const double alpha, double *x, const double beta, double *y, const int nrhs)
This module contains the generic interfaces for all the analysis routines.
void dqrm_residual_orth_c(struct dqrm_spmat_type_c *qrm_spmat_c, double *r, const int nrhs, double *nrm)
void dqrm_spmat_destroy_c(struct dqrm_spmat_type_c *qrm_spmat_c)
subroutine dqrm_matnrm(qrm_mat, ntype, nrm)
This subroutine computes the matrix norm. The return value is a real scalar.
subroutine dqrm_factorize(qrm_mat, transp)
This routine is the main factorization driver.