QR_MUMPS
qrm_spmat_mod.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 
37 
41  use qrm_common_mod
42  use qrm_error_mod
43  use qrm_mem_mod
44  use qrm_adata_mod
45  use _qrm_fdata_mod
46  use qrm_error_mod
47 
49  interface qrm_spmat_alloc
50  module procedure _qrm_spmat_alloc
51  end interface
52 
54  interface qrm_spmat_init
55  module procedure _qrm_spmat_init
56  end interface
57 
59  interface qrm_cntl_init
60  module procedure _qrm_cntl_init
61  end interface
62 
65  module procedure _qrm_spmat_convert
66  end interface
67 
69  interface qrm_spmat_copy
70  module procedure _qrm_spmat_copy
71  end interface
72 
75  module procedure _qrm_spmat_destroy
76  end interface
77 
78 
82  ! ! !! @link ::_qrm_psetl @endlink routines
83  !!
84  interface qrm_set
85  module procedure _qrm_pseti, _qrm_psetr!, _qrm_psetl
86  end interface
87 
91  ! !! @link ::_qrm_pgetl @endlink routines
92  !!
93  interface qrm_get
94  module procedure _qrm_pgeti, _qrm_pgetii, _qrm_pgetr!, _qrm_pgetl
95  end interface
96 
97 
99  interface qrm_get_r
100  module procedure _qrm_get_r
101  end interface
102 
104 
142 
143  integer :: icntl(20)=0
154  real(kind(1.d0)) :: rcntl(10)=0.d0
162  integer(kind=8) :: gstats(10)=0
164  integer, pointer, dimension(:) :: iptr => null()
166  integer, pointer, dimension(:) :: jptr => null()
168  integer, pointer, dimension(:) :: irn => null()
170  integer, pointer, dimension(:) :: jcn => null()
172  _qrm_data, pointer, dimension(:) :: val => null()
174  integer :: m=0
176  integer :: n=0
178  integer :: nz=0
181  integer, pointer, dimension(:) :: cperm_in => null()
184  type(qrm_adata_type) :: adata
187  type(_qrm_fdata_type) :: fdata
189  character(len=3) :: fmt='coo'
190  end type _qrm_spmat_type
191 
192 
193 contains
194 
216  subroutine _qrm_spmat_alloc(qrm_spmat, nz, m, n, fmt)
218  use qrm_error_mod
219  implicit none
220 
221  type(_qrm_spmat_type), intent(inout) :: qrm_spmat
222  integer, intent(in) :: nz, m, n
223  character, intent(in) :: fmt*(*)
224 
225  logical :: lsamen ! LAPACK subroutine to test strings
226  ! error management
227  integer :: err_act
228  character(len=*), parameter :: name='_qrm_spmat_alloc'
229 
230  call qrm_err_act_save(err_act)
231 
232 #if defined(debug)
233  __qrm_prnt_dbg('("Allocating Matrix")')
234 #endif
235 
236  if(fmt .eq. 'coo') then
237  call qrm_palloc(qrm_spmat%irn, nz)
238  call qrm_palloc(qrm_spmat%jcn, nz)
239  call qrm_palloc(qrm_spmat%val, nz)
240  __qrm_check_ret(name,'qrm_palloc',9999)
241  else if(fmt .eq. 'csr') then
242  call qrm_palloc(qrm_spmat%iptr, m+1)
243  call qrm_palloc(qrm_spmat%jcn , nz)
244  call qrm_palloc(qrm_spmat%val , nz)
245  __qrm_check_ret(name,'qrm_palloc',9999)
246  else if(fmt .eq. 'csc') then
247  call qrm_palloc(qrm_spmat%irn , nz)
248  call qrm_palloc(qrm_spmat%jptr, n+1)
249  call qrm_palloc(qrm_spmat%val , nz)
250  __qrm_check_ret(name,'qrm_palloc',9999)
251  else
252  call qrm_err_push(1, '_qrm_spmat_convert',aed=fmt)
253  goto 9999
254  end if
255 
256  qrm_spmat%nz = nz
257  qrm_spmat%m = m
258  qrm_spmat%n = n
259 
260  call qrm_err_act_restore(err_act)
261  return
262 
263 9999 continue ! error management
264  call qrm_err_act_restore(err_act)
265  if(err_act .eq. qrm_abort_) then
266  call qrm_err_check()
267  end if
268  return
269 
270  end subroutine _qrm_spmat_alloc
271 
272 
278  subroutine _qrm_spmat_init(qrm_spmat)
280  implicit none
281  type(_qrm_spmat_type), intent(inout) :: qrm_spmat
282 
283  character(LEN=10) :: str
284  integer :: ierr
285 
286  call _qrm_cntl_init(qrm_spmat)
287 
288 
289  nullify(qrm_spmat%iptr, qrm_spmat%jptr, qrm_spmat%irn, qrm_spmat%jcn, &
290  & qrm_spmat%val, qrm_spmat%cperm_in)
291 
292  return
293 
294  end subroutine _qrm_spmat_init
295 
296 
302  subroutine _qrm_cntl_init(qrm_spmat)
304  use qrm_common_mod
305  implicit none
306  type(_qrm_spmat_type), intent(inout) :: qrm_spmat
307 
308  character(LEN=10) :: str
309  integer :: ierr
310 
311  ! set default values for icntl and rcntl
312  qrm_spmat%icntl(qrm_ordering_) = qrm_auto_
313  qrm_spmat%icntl(qrm_minamalg_) = 4
314  qrm_spmat%icntl(qrm_nb_) = 120
315  qrm_spmat%icntl(qrm_keeph_) = qrm_yes_
316  qrm_spmat%icntl(qrm_ib_) = 120
317  qrm_spmat%icntl(qrm_rhsnb_) = -1
318  qrm_spmat%icntl(qrm_rhsnthreads_) = 1
319  qrm_spmat%icntl(qrm_nlz_) = 8
320  qrm_spmat%icntl(qrm_cnode_) = 1
321  qrm_spmat%icntl(qrm_sing_) = qrm_no_
322 
323  qrm_spmat%rcntl(qrm_amalgth_) = 0.05
324  qrm_spmat%rcntl(qrm_rweight_) = 0.001
325  qrm_spmat%fmt = 'coo'
326 
327  call get_environment_variable(name="QRM_NUM_THREADS",value=str, status=ierr)
328  if(ierr .eq. 1) then
329  qrm_spmat%icntl(qrm_nthreads_) = 1
330  else
331  read(str,*)qrm_spmat%icntl(qrm_nthreads_)
332  end if
333 
334  return
335 
336  end subroutine _qrm_cntl_init
337 
338 
339 
340 
341 
355  subroutine _qrm_spmat_convert(in_mat, out_mat, fmt, values)
356  implicit none
357 
358  type(_qrm_spmat_type), intent(in) :: in_mat
359  type(_qrm_spmat_type) :: out_mat
360  character, intent(in) :: fmt*(*)
361  logical, optional :: values
362  ! error management
363  integer :: err_act
364  character(len=*), parameter :: name='_qrm_spmat_convert'
365 
366  call qrm_err_act_save(err_act)
367 
368  select case(in_mat%fmt)
369  case('csc')
370  select case(fmt)
371  case('csr')
372  call _qrm_csc_to_csr(in_mat, out_mat, values)
373  case default
374  call qrm_err_push(1, name,aed=in_mat%fmt)
375  goto 9999
376  end select
377  case('coo')
378  select case(fmt)
379  case('csc')
380  call _qrm_coo_to_csc(in_mat, out_mat, values)
381  case default
382  call qrm_err_push(1, name,aed=in_mat%fmt)
383  goto 9999
384  end select
385  case default
386  call qrm_err_push(1, name,aed=in_mat%fmt)
387  goto 9999
388  end select
389 
390  out_mat%icntl = in_mat%icntl
391  out_mat%rcntl = in_mat%rcntl
392 
393  call qrm_err_act_restore(err_act)
394  return
395 
396 9999 continue ! error management
397  call qrm_err_act_restore(err_act)
398  if(err_act .eq. qrm_abort_) then
399  call qrm_err_check()
400  end if
401  return
402 
403  end subroutine _qrm_spmat_convert
404 
405 
417  subroutine _qrm_coo_to_csc(in_mat, out_mat, values)
419  implicit none
420 
421  type(_qrm_spmat_type), intent(in) :: in_mat
422  type(_qrm_spmat_type) :: out_mat
423  logical, optional :: values
424 
425  integer, allocatable :: work(:)
426  logical :: ivalues, ob
427  integer :: i, j, idx, k, m, n
428  ! error management
429  integer :: err_act
430  character(len=*), parameter :: name='_qrm_coo_to_csc'
431 
432  call qrm_err_act_save(err_act)
433 
434  if(present(values)) then
435  ivalues = values
436  else
437  ivalues = .true.
438  end if
439 
440  call qrm_aalloc(work, in_mat%n+1)
441 
442  call qrm_prealloc(out_mat%jptr, in_mat%n+1)
443  call qrm_prealloc(out_mat%irn , in_mat%nz)
444  if(ivalues) call qrm_prealloc(out_mat%val , in_mat%nz)
445  __qrm_check_ret(name,'qrm_alloc',9999)
446 
447  work=0
448  ob = .false.
449 
450  m = in_mat%m
451  n = in_mat%n
452 
453  ! first loop to calculate # of nz per column
454  do k=1, in_mat%nz
455  j = in_mat%jcn(k)
456  i = in_mat%irn(k)
457  if((j.gt.0) .and. (j.le. n) .and. (i.gt.0) .and. (i.le. m) ) then
458  work(j) = work(j)+1
459  else
460  ! out of bounds coefficients. ignore and print a warning at the end
461  ob = .true.
462  end if
463  end do
464 
465  if(ob) then
466  __qrm_prnt_dbg('("** Out-of-bounds coefficients present **")')
467  end if
468 
469  ! loop to convert the counts into ptrs
470  out_mat%jptr(1) = 1
471  do j=2, n+1
472  out_mat%jptr(j) = out_mat%jptr(j-1)+work(j-1)
473  end do
474 
475 
476  ! last loop to put things into place
477  work=0
478  ! instead of putting an "if" inside the loop
479  ! I put it here to gain some speed
480  if(ivalues) then
481  do k=1, in_mat%nz
482  j = in_mat%jcn(k)
483  i = in_mat%irn(k)
484  if((j.le.0) .or. (j.gt. n) .or. (i.le.0) .or. (i.gt. m) ) cycle
485  idx = out_mat%jptr(j)+work(j)
486  out_mat%irn(idx) = i
487  out_mat%val(idx) = in_mat%val(k)
488  work(j) = work(j)+1
489  end do
490  else
491  do k=1, in_mat%nz
492  j = in_mat%jcn(k)
493  i = in_mat%irn(k)
494  if((j.le.0) .or. (j.gt. n) .or. (i.le.0) .or. (i.gt. m) ) cycle
495  idx = out_mat%jptr(j)+work(j)
496  out_mat%irn(idx) = i
497  work(j) = work(j)+1
498  end do
499  end if
500 
501  call qrm_adealloc(work)
502  __qrm_check_ret(name,'qrm_adelloc',9999)
503 
504  out_mat%m = in_mat%m
505  out_mat%n = in_mat%n
506  out_mat%nz = in_mat%nz
507  out_mat%fmt = 'csc'
508 
509  call qrm_err_act_restore(err_act)
510  return
511 
512 9999 continue ! error management
513  call qrm_err_act_restore(err_act)
514  if(err_act .eq. qrm_abort_) then
515  call qrm_err_check()
516  end if
517  return
518 
519  end subroutine _qrm_coo_to_csc
520 
532  subroutine _qrm_csc_to_csr(in_mat, out_mat, values)
534  type(_qrm_spmat_type), intent(in) :: in_mat
535  type(_qrm_spmat_type) :: out_mat
536  logical, optional :: values
537 
538  integer, allocatable :: work(:)
539 
540  logical :: ivalues, ob
541  integer :: i, j, idx, ii, m, n
542  ! error management
543  integer :: err_act
544  character(len=*), parameter :: name='_qrm_csc_to_csr'
545 
546  call qrm_err_act_save(err_act)
547 
548  if(present(values)) then
549  ivalues=values
550  else
551  ivalues = .true.
552  end if
553 
554  ob = .false.
555 
556  m = in_mat%m
557  n = in_mat%n
558 
559  call qrm_aalloc(work, m+1)
560 
561  call qrm_prealloc(out_mat%iptr, m+1)
562  call qrm_prealloc(out_mat%jcn , in_mat%nz)
563  if(ivalues) call qrm_prealloc(out_mat%val , in_mat%nz)
564  __qrm_check_ret(name,'qrm_alloc',9999)
565 
566  work=0
567  ! first loop to calculate # of nz per row
568  do j = 1, n
569  do ii= in_mat%jptr(j), in_mat%jptr(j+1)-1
570  i = in_mat%irn(ii)
571  if((i.gt.0) .and. (i.le.m)) then
572  work(i) = work(i)+1
573  else
574  ob = .true.
575  end if
576  end do
577  end do
578 
579  if(ob) then
580  __qrm_prnt_dbg('("** Out-of-bounds coefficients present **")')
581  end if
582 
583  ! loop to convert the counts into ptrs
584  out_mat%iptr(1) = 1
585  do j=2, m+1
586  out_mat%iptr(j) = out_mat%iptr(j-1)+work(j-1)
587  end do
588 
589 
590  ! last loop to put things into place
591  work=0
592  ! instead of putting an "if" inside the loop
593  ! I put it here to gain some speed
594  if(ivalues) then
595  do j = 1, n
596  do ii= in_mat%jptr(j), in_mat%jptr(j+1)-1
597  i = in_mat%irn(ii)
598  if((i.le.0) .or. (i.gt.m)) cycle
599  idx = out_mat%iptr(i)+work(i)
600  out_mat%jcn(idx) = j
601  out_mat%val(idx) = in_mat%val(ii)
602  work(i) = work(i)+1
603  end do
604  end do
605  else
606  do j = 1, n
607  do ii= in_mat%jptr(j), in_mat%jptr(j+1)-1
608  i = in_mat%irn(ii)
609  if((i.le.0) .or. (i.gt.m)) cycle
610  idx = out_mat%iptr(i)+work(i)
611  out_mat%jcn(idx) = j
612  work(i) = work(i)+1
613  end do
614  end do
615  end if
616 
617  call qrm_adealloc(work)
618  __qrm_check_ret(name,'qrm_adelloc',9999)
619 
620  out_mat%m = in_mat%m
621  out_mat%n = in_mat%n
622  out_mat%nz = in_mat%nz
623  out_mat%fmt = 'csr'
624 
625  call qrm_err_act_restore(err_act)
626  return
627 
628 9999 continue ! error management
629  call qrm_err_act_restore(err_act)
630  if(err_act .eq. qrm_abort_) then
631  call qrm_err_check()
632  end if
633  return
634 
635  end subroutine _qrm_csc_to_csr
636 
637 
638 
639 
640 
641 
652  subroutine _qrm_spmat_copy(in_mat, out_mat, values)
654  type(_qrm_spmat_type), intent(in) :: in_mat
655  type(_qrm_spmat_type) :: out_mat
656  logical, optional :: values
657 
658  logical :: ivalues=.true.
659  ! error management
660  integer :: err_act
661  character(len=*), parameter :: name='_qrm_spmat_copy'
662 
663  call qrm_err_act_save(err_act)
664 
665  ! TODO complete with other types
666 
667  if(present(values)) ivalues=values
668 
669  select case(in_mat%fmt)
670  case('csc')
671  call qrm_prealloc(out_mat%jptr, in_mat%n+1)
672  call qrm_prealloc(out_mat%irn, in_mat%nz)
673  __qrm_check_ret(name,'qrm_prelloc',9999)
674 
675  do i=1, in_mat%n+1
676  out_mat%jptr(i) = in_mat%jptr(i)
677  end do
678  do i=1, in_mat%nz
679  out_mat%irn(i) = in_mat%irn(i)
680  end do
681  if(ivalues) then
682  call qrm_prealloc(out_mat%val, in_mat%nz)
683  __qrm_check_ret(name,'qrm_prealloc',9999)
684  out_mat%val = in_mat%val
685  end if
686  case('coo')
687  call qrm_prealloc(out_mat%jcn, in_mat%nz)
688  call qrm_prealloc(out_mat%irn, in_mat%nz)
689  __qrm_check_ret(name,'qrm_prealloc',9999)
690  do i=1, in_mat%nz
691  out_mat%jcn(i) = in_mat%jcn(i)
692  out_mat%irn(i) = in_mat%irn(i)
693  end do
694  if(ivalues) then
695  call qrm_prealloc(out_mat%val, in_mat%nz)
696  __qrm_check_ret(name,'qrm_realloc',9999)
697  out_mat%val = in_mat%val
698  end if
699  case default
700  call qrm_err_push(1, name,aed=in_mat%fmt)
701  goto 9999
702  end select
703 
704  out_mat%n = in_mat%n
705  out_mat%m = in_mat%m
706  out_mat%nz = in_mat%nz
707  out_mat%fmt = in_mat%fmt
708  out_mat%icntl = in_mat%icntl
709  out_mat%rcntl = in_mat%rcntl
710 
711  call qrm_err_act_restore(err_act)
712  return
713 
714 9999 continue ! error management
715  call qrm_err_act_restore(err_act)
716  if(err_act .eq. qrm_abort_) then
717  call qrm_err_check()
718  end if
719  return
720 
721  end subroutine _qrm_spmat_copy
722 
729  subroutine _qrm_spmat_destroy(qrm_spmat, all)
731  use qrm_mem_mod
732  implicit none
733 
734  type(_qrm_spmat_type) :: qrm_spmat
735  logical, optional :: all
736 
737  logical :: iall
738  ! error management
739  integer :: err_act
740  character(len=*), parameter :: name='_qrm_spmat_destroy'
741 
742  call qrm_err_act_save(err_act)
743 
744  if(present(all)) then
745  iall = all
746  else
747  iall = .false.
748  end if
749 
750  if(iall) then
751  call qrm_pdealloc(qrm_spmat%irn)
752  call qrm_pdealloc(qrm_spmat%jcn)
753  call qrm_pdealloc(qrm_spmat%iptr)
754  call qrm_pdealloc(qrm_spmat%jptr)
755  call qrm_pdealloc(qrm_spmat%val)
756  call qrm_pdealloc(qrm_spmat%cperm_in)
757  __qrm_check_ret(name,'qrm_pdealloc',9999)
758  end if
759 
760  qrm_spmat%n = 0
761  qrm_spmat%m = 0
762  qrm_spmat%nz = 0
763  qrm_spmat%fmt = ''
764 
765  call qrm_adata_destroy(qrm_spmat%adata)
766  __qrm_check_ret(name,name,9999)
767  call _qrm_fdata_destroy(qrm_spmat%fdata)
768  __qrm_check_ret(name,name,9999)
769 
770  call qrm_err_act_restore(err_act)
771  return
772 
773 9999 continue ! error management
774  call qrm_err_act_restore(err_act)
775  if(err_act .eq. qrm_abort_) then
776  call qrm_err_check()
777  end if
778  return
779 
780  end subroutine _qrm_spmat_destroy
781 
782 
783 
784  ! The following subroutine set or get control parameters from the
785  ! cntl or rcntl control arrays. All the set and get routines are
786  ! gathered under the same, overloaded interface, respectively
787 
788 
823  subroutine _qrm_pseti(qrm_spmat, string, ival)
825  use qrm_string_mod
826  use qrm_error_mod
827  implicit none
828 
829  type(_qrm_spmat_type) :: qrm_spmat
830  character(len=*) :: string
831  integer :: ival
832 
833  character(len=len(string)) :: istring
834  ! error management
835  integer :: err_act
836  character(len=*), parameter :: name='_qrm_pseti'
837 
838  call qrm_err_act_save(err_act)
839 
840  istring = qrm_str_tolower(string)
841  if(index(istring,'qrm_ordering') .eq. 1) then
842  qrm_spmat%icntl(qrm_ordering_) = ival
843  else if (index(istring,'qrm_minamalg') .eq. 1) then
844  qrm_spmat%icntl(qrm_minamalg_) = ival
845  else if (index(istring,'qrm_nb') .eq. 1) then
846  qrm_spmat%icntl(qrm_nb_) = ival
847  if (qrm_spmat%icntl(qrm_ib_).gt.qrm_spmat%icntl(qrm_nb_)) then
848  __qrm_prnt_msg('("Warning: qrm_ib is being set equal to qrm_nb")')
849  qrm_spmat%icntl(qrm_ib_) = qrm_spmat%icntl(qrm_nb_)
850  end if
851  else if (index(istring,'qrm_ib') .eq. 1) then
852  qrm_spmat%icntl(qrm_ib_) = ival
853  if (qrm_spmat%icntl(qrm_ib_).gt.qrm_spmat%icntl(qrm_nb_)) then
854  __qrm_prnt_msg('("Warning: qrm_nb is being set equal to qrm_ib")')
855  qrm_spmat%icntl(qrm_nb_) = qrm_spmat%icntl(qrm_ib_)
856  end if
857  else if (index(istring,'qrm_rhsnb') .eq. 1) then
858  qrm_spmat%icntl(qrm_rhsnb_) = ival
859  else if (index(istring,'qrm_nthreads') .eq. 1) then
860  qrm_spmat%icntl(qrm_nthreads_) = ival
861  else if (index(istring,'qrm_rhsnthreads') .eq. 1) then
862  qrm_spmat%icntl(qrm_rhsnthreads_) = ival
863  else if (index(istring,'qrm_keeph') .eq. 1) then
864  if(ival .eq. qrm_yes_) then
865  qrm_spmat%icntl(qrm_keeph_) = ival
866  else
867  qrm_spmat%icntl(qrm_keeph_) = qrm_no_
868  end if
869  else if (index(istring,'qrm_sing') .eq. 1) then
870  if(ival .eq. qrm_yes_) then
871  qrm_spmat%icntl(qrm_sing_) = ival
872  else
873  qrm_spmat%icntl(qrm_sing_) = qrm_no_
874  end if
875  else if (index(istring,'qrm_nlz') .eq. 1) then
876  qrm_spmat%icntl(qrm_nlz_) = ival
877  else if (index(istring,'qrm_cnode') .eq. 1) then
878  qrm_spmat%icntl(qrm_cnode_) = ival
879  else
880  call qrm_err_push(23, name, aed=string)
881  goto 9999
882  end if
883 
884  call qrm_err_act_restore(err_act)
885  return
886 
887 9999 continue ! error management
888  call qrm_err_act_restore(err_act)
889  if(err_act .eq. qrm_abort_) then
890  call qrm_err_check()
891  end if
892 
893  return
894 
895  end subroutine _qrm_pseti
896 
897 
910  subroutine _qrm_psetr(qrm_spmat, string, rval)
912  use qrm_common_mod
913  use qrm_string_mod
914  use qrm_error_mod
915  implicit none
916 
917  type(_qrm_spmat_type) :: qrm_spmat
918  character(len=*) :: string
919  real(kind(1.d0)) :: rval
920 
921  character(len=len(string)) :: istring
922  ! error management
923  integer :: err_act
924  character(len=*), parameter :: name='_qrm_psetr'
925 
926  call qrm_err_act_save(err_act)
927 
928  istring = qrm_str_tolower(string)
929 
930  if(index(istring,'qrm_amalgth') .eq. 1) then
931  qrm_spmat%rcntl(qrm_amalgth_) = rval
932  else if(index(istring,'qrm_rweight') .eq. 1) then
933  qrm_spmat%rcntl(qrm_rweight_) = rval
934  else
935  call qrm_err_push(23, name, aed=string)
936  goto 9999
937  end if
938 
939  call qrm_err_act_restore(err_act)
940  return
941 
942 9999 continue ! error management
943  call qrm_err_act_restore(err_act)
944  if(err_act .eq. qrm_abort_) then
945  call qrm_err_check()
946  end if
947 
948  return
949 
950  end subroutine _qrm_psetr
951 
952 
953 
954  ! subroutine _qrm_psetl(qrm_spmat, string, lval)
955 
956  ! use qrm_string_mod
957  ! use qrm_error_mod
958  ! implicit none
959 
960  ! type(_qrm_spmat_type) :: qrm_spmat
961  ! character(len=*) :: string
962  ! logical :: lval
963 
964  ! character(len=len(string)) :: istring
965  ! ! error management
966  ! integer :: err_act
967  ! character(len=*), parameter :: name='_qrm_psetl'
968 
969  ! call qrm_err_act_save(err_act)
970 
971  ! istring = qrm_str_tolower(string)
972 
973  ! if(index(istring,'qrm_keeph') .eq. qrm_yes_) then
974  ! if(lval) then
975  ! qrm_spmat%icntl(qrm_keeph_) = qrm_yes_
976  ! else
977  ! qrm_spmat%icntl(qrm_keeph_) = qrm_no_
978  ! end if
979  ! else if(index(istring,'qrm_sing') .eq. 1) then
980  ! if(lval) then
981  ! qrm_spmat%icntl(qrm_sing_) = 1
982  ! else
983  ! qrm_spmat%icntl(qrm_sing_) = 0
984  ! end if
985  ! else
986  ! call qrm_err_push(23, name, aed=string)
987  ! goto 9999
988  ! end if
989 
990  ! call qrm_err_act_restore(err_act)
991  ! return
992 
993 ! 9999 continue ! error management
994  ! call qrm_err_act_restore(err_act)
995  ! if(err_act .eq. qrm_abort_) then
996  ! call qrm_err_check()
997  ! end if
998 
999  ! return
1000 
1001  ! end subroutine _qrm_psetl_
1002 
1003 
1007  subroutine _qrm_pgeti(qrm_spmat, string, ival)
1009  use qrm_string_mod
1010  use qrm_error_mod
1011  implicit none
1012 
1013  type(_qrm_spmat_type) :: qrm_spmat
1014  character(len=*) :: string
1015  integer :: ival
1016 
1017  character(len=len(string)) :: istring
1018  integer(kind=8) :: iival
1019  ! error management
1020  integer :: err_act
1021  character(len=*), parameter :: name='_qrm_pgeti'
1022 
1023  call qrm_err_act_save(err_act)
1024 
1025  call _qrm_pgetii(qrm_spmat, string, iival)
1026  __qrm_check_ret(name,'qrm_pgetii',9999)
1027 
1028  ival = iival
1029 
1030  call qrm_err_act_restore(err_act)
1031  return
1032 
1033 9999 continue ! error management
1034  call qrm_err_act_restore(err_act)
1035  if(err_act .eq. qrm_abort_) then
1036  call qrm_err_check()
1037  end if
1038 
1039  return
1040 
1041  end subroutine _qrm_pgeti
1042 
1046  subroutine _qrm_pgetii(qrm_spmat, string, ival)
1048  use qrm_string_mod
1049  use qrm_error_mod
1050  implicit none
1051 
1052  type(_qrm_spmat_type) :: qrm_spmat
1053  character(len=* ) :: string
1054  integer(kind=8) :: ival
1055 
1056  character(len=len(string)) :: istring
1057  ! error management
1058  integer :: err_act
1059  character(len=*), parameter :: name='_qrm_pgetii'
1060 
1061  call qrm_err_act_save(err_act)
1062 
1063  istring = qrm_str_tolower(string)
1064 
1065  if(index(istring,'qrm_ordering') .eq. 1) then
1066  ival = qrm_spmat%icntl(qrm_ordering_)
1067  else if (index(istring,'qrm_minamalg') .eq. 1) then
1068  ival = qrm_spmat%icntl(qrm_minamalg_)
1069  else if (index(istring,'qrm_nb') .eq. 1) then
1070  ival = qrm_spmat%icntl(qrm_nb_)
1071  else if (index(istring,'qrm_ib') .eq. 1) then
1072  ival = qrm_spmat%icntl(qrm_ib_)
1073  else if (index(istring,'qrm_rhsnb') .eq. 1) then
1074  ival = qrm_spmat%icntl(qrm_rhsnb_)
1075  else if (index(istring,'qrm_nthreads') .eq. 1) then
1076  ival = qrm_spmat%icntl(qrm_nthreads_)
1077  else if (index(istring,'qrm_rhsnthreads') .eq. 1) then
1078  ival = qrm_spmat%icntl(qrm_rhsnthreads_)
1079  else if (index(istring,'qrm_keeph') .eq. 1) then
1080  ival = qrm_spmat%icntl(qrm_keeph_)
1081  else if (index(istring,'qrm_sing') .eq. 1) then
1082  ival = qrm_spmat%icntl(qrm_sing_)
1083  else if (index(istring,'qrm_e_nnz_r') .eq. 1) then
1084  ival = qrm_spmat%gstats(qrm_e_nnz_r_)
1085  else if (index(istring,'qrm_e_nnz_h') .eq. 1) then
1086  ival = qrm_spmat%gstats(qrm_e_nnz_h_)
1087  else if (index(istring,'qrm_e_facto_flops') .eq. 1) then
1088  ival = qrm_spmat%gstats(qrm_e_facto_flops_)
1089  else if (index(istring,'qrm_nnz_r') .eq. 1) then
1090  ival = qrm_spmat%gstats(qrm_nnz_r_)
1091  else if (index(istring,'qrm_nnz_h') .eq. 1) then
1092  ival = qrm_spmat%gstats(qrm_nnz_h_)
1093  else if (index(istring,'qrm_facto_flops') .eq. 1) then
1094  ival = qrm_spmat%gstats(qrm_facto_flops_)
1095  else
1096  call qrm_err_push(23, name, aed=string)
1097  goto 9999
1098  end if
1099 
1100  call qrm_err_act_restore(err_act)
1101  return
1102 
1103 9999 continue ! error management
1104  call qrm_err_act_restore(err_act)
1105  if(err_act .eq. qrm_abort_) then
1106  call qrm_err_check()
1107  end if
1108 
1109  return
1110 
1111  end subroutine _qrm_pgetii
1112 
1113 
1114 
1119  subroutine _qrm_pgetr(qrm_spmat, string, rval)
1121  use qrm_common_mod
1122  use qrm_string_mod
1123  use qrm_error_mod
1124  implicit none
1125 
1126  type(_qrm_spmat_type) :: qrm_spmat
1127  character(len=*) :: string
1128  real(kind(1.d0)) :: rval
1129 
1130  character(len=len(string)) :: istring
1131  ! error management
1132  integer :: err_act
1133  character(len=*), parameter :: name='_qrm_pgetr'
1134 
1135  call qrm_err_act_save(err_act)
1136 
1137  istring = qrm_str_tolower(string)
1138 
1139  if(index(istring,'qrm_amalgth') .eq. 1) then
1140  rval = qrm_spmat%rcntl(qrm_amalgth_)
1141  else
1142  call qrm_err_push(23, name, aed=string)
1143  goto 9999
1144  end if
1145 
1146  call qrm_err_act_restore(err_act)
1147  return
1148 
1149 9999 continue ! error management
1150  call qrm_err_act_restore(err_act)
1151  if(err_act .eq. qrm_abort_) then
1152  call qrm_err_check()
1153  end if
1154 
1155  return
1156 
1157  end subroutine _qrm_pgetr
1158 
1159 
1160  ! subroutine _qrm_pgetl(qrm_spmat, string, lval)
1161 
1162  ! use qrm_string_mod
1163  ! use qrm_error_mod
1164  ! implicit none
1165 
1166  ! type(_qrm_spmat_type) :: qrm_spmat
1167  ! character(len=*) :: string
1168  ! logical :: lval
1169 
1170  ! character(len=len(string)) :: istring
1171  ! ! error management
1172  ! integer :: err_act
1173  ! character(len=*), parameter :: name='_qrm_pgetl'
1174 
1175  ! call qrm_err_act_save(err_act)
1176 
1177  ! istring = qrm_str_tolower(string)
1178 
1179  ! if(index(istring,'qrm_keeph') .eq. 1) then
1180  ! if (qrm_spmat%icntl(qrm_keeph_) .eq. qrm_yes_) then
1181  ! lval = .true.
1182  ! else
1183  ! lval = .false.
1184  ! end if
1185  ! else if(index(istring,'qrm_sing') .eq. 1) then
1186  ! if(qrm_spmat%icntl(qrm_sing_) .eq. 1) then
1187  ! lval = .true.
1188  ! else
1189  ! lval = .false.
1190  ! end if
1191  ! else
1192  ! call qrm_err_push(23, name, aed=string)
1193  ! goto 9999
1194  ! end if
1195 
1196  ! call qrm_err_act_restore(err_act)
1197  ! return
1198 
1199 ! 9999 continue ! error management
1200  ! call qrm_err_act_restore(err_act)
1201  ! if(err_act .eq. qrm_abort_) then
1202  ! call qrm_err_check()
1203  ! end if
1204 
1205  ! return
1206 
1207  ! end subroutine _qrm_pgetl
1208 
1209 
1212  subroutine _qrm_check_spmat(qrm_spmat, op)
1214  use qrm_string_mod
1215  use qrm_error_mod
1216  implicit none
1217 
1218  type(_qrm_spmat_type) :: qrm_spmat
1219  integer, optional :: op
1220 
1221  integer :: iop
1222 
1223  ! error management
1224  integer :: err_act
1225  character(len=*), parameter :: name='_qrm_check_spmat'
1226 
1227  call qrm_err_act_save(err_act)
1228 
1229  if(present(op)) then
1230  iop = op
1231  else
1232  iop = qrm_allop_
1233  end if
1234 
1235  if((qrm_spmat%m .lt. 0) .or. (qrm_spmat%n .lt. 0) .or. &
1236  & (qrm_spmat%nz .lt. 0) .or. &
1237  & (qrm_spmat%nz .gt. (int(qrm_spmat%n,kind=8)*int(qrm_spmat%m,kind=8)))) then
1238  call qrm_err_push(29, name,ied=(/qrm_spmat%m,qrm_spmat%n,qrm_spmat%nz,0,0/))
1239  goto 9999
1240  end if
1241 
1242 
1243  if((iop.eq.qrm_allop_) .or. (iop.eq.qrm_analyse_)) then
1244 
1245  ! all the potential cases of conflict with the orderings
1246  select case(qrm_spmat%icntl(qrm_ordering_))
1247  case(:-1,6:)
1248  call qrm_err_push(9, name,ied=(/qrm_spmat%icntl(qrm_ordering_),0,0,0,0/))
1249  goto 9999
1250  case (qrm_given_)
1251  if(qrm_spmat%icntl(qrm_sing_) .eq. qrm_yes_) then
1252  call qrm_err_push(27, name,ied=(/qrm_ordering_,qrm_sing_,0,0,0/))
1253  goto 9999
1254  end if
1255  end select
1256 
1257  ! all the potential cases of conflict with the orderings
1258  select case(qrm_spmat%icntl(qrm_nb_))
1259  case(:-1)
1260  call qrm_err_push(28, name, ied=(/qrm_spmat%icntl(qrm_nb_),0,0,0,0/))
1261  goto 9999
1262  case default
1263  if(qrm_spmat%icntl(qrm_nb_) .lt. qrm_spmat%icntl(qrm_ib_)) then
1264  call qrm_err_push(27, name,ied=(/qrm_nb_,qrm_ib_,0,0,0/))
1265  goto 9999
1266  end if
1267  end select
1268 
1269  select case(qrm_spmat%icntl(qrm_ib_))
1270  case(:-1)
1271  call qrm_err_push(28, name, ied=(/qrm_spmat%icntl(qrm_ib_),0,0,0,0/))
1272  goto 9999
1273  end select
1274  end if
1275 
1276  call qrm_err_act_restore(err_act)
1277  return
1278 
1279 9999 continue ! error management
1280  call qrm_err_act_restore(err_act)
1281  if(err_act .eq. qrm_abort_) then
1282  call qrm_err_check()
1283  end if
1284  return
1285 
1286  end subroutine _qrm_check_spmat
1287 
1288 
1289  subroutine _qrm_get_r(qrm_mat, r)
1291  use qrm_error_mod
1292  use _qrm_fdata_mod
1293  use qrm_mem_mod
1294  implicit none
1295 
1296  type(_qrm_spmat_type), target :: qrm_mat
1297  type(_qrm_spmat_type) :: r
1298  type(_qrm_front_type), pointer :: front
1299 
1300  integer :: cnt, fcnt, f, jp, pk, j, k, n, c, rbcnt, rtcnt, i, rps, du, eu
1301 
1302  r%nz = qrm_mat%gstats(qrm_nnz_r_)
1303  r%m = size(qrm_mat%adata%rperm)
1304  r%n = size(qrm_mat%adata%cperm)
1305 
1306  call qrm_palloc(r%irn, r%nz)
1307  call qrm_palloc(r%jcn, r%nz)
1308  call qrm_palloc(r%val, r%nz)
1309  call qrm_aalloc(r%adata%rperm, r%m)
1310  call qrm_aalloc(r%adata%cperm, r%n)
1311  r%adata%cperm = qrm_mat%adata%cperm
1312 
1313  rtcnt = 1
1314  rbcnt = min(r%m,r%n)+1
1315  rps = 0
1316  cnt = 1
1317  do f = 1, qrm_mat%adata%nnodes
1318  front => qrm_mat%fdata%front_list(f)
1319  rps = rps + front%npiv + front%m-front%ne
1320  r%adata%rperm(rtcnt:rtcnt+front%npiv-1) = front%rows(1:front%npiv)
1321  r%adata%rperm(rbcnt:rbcnt + front%m-front%ne-1) = front%rows(front%ne+1:front%m)
1322  rtcnt = rtcnt+front%npiv
1323  rbcnt = rbcnt + front%m-front%ne
1324 
1325  fcnt = 1
1326  outer: do jp = 1, front%npiv, front%nb
1327  pk = min(front%nb, front%npiv-jp+1)
1328  if(pk .le. 0) exit outer
1329 
1330  inner: do j = jp, jp+pk-1, front%ib
1331  k = min(front%ib, jp+pk - j)
1332  if(k .le. 0) exit inner
1333  n = front%n-j+1
1334 
1335  do c=1, k
1336  r%irn(cnt:cnt+c-1) = front%rows(j:j+c-1)
1337  r%jcn(cnt:cnt+c-1) = front%cols(j+c-1)
1338  r%val(cnt:cnt+c-1) = front%r(fcnt:fcnt+c-1)
1339  fcnt = fcnt+k
1340  cnt = cnt+c
1341  end do
1342 
1343  do c=k+1, n
1344  r%irn(cnt:cnt+k-1) = front%rows(j:j+k-1)
1345  r%jcn(cnt:cnt+k-1) = front%cols(j+c-1)
1346  r%val(cnt:cnt+k-1) = front%r(fcnt:fcnt+k-1)
1347  fcnt = fcnt+k
1348  cnt = cnt+k
1349  end do
1350 
1351  end do inner
1352  end do outer
1353 
1354  end do
1355 
1356  if(rbcnt .ne. r%m+1) then
1357  __qrm_prnt_dbg('("_qrm_get_r -- The matrix contains empty rows")')
1358  r%adata%rperm(rbcnt:r%m) = qrm_mat%adata%rperm(rbcnt:r%m)
1359  end if
1360 
1361  if(rtcnt.lt.min(r%m,r%n)) then
1362  __qrm_prnt_err('("_qrm_get_r -- The R matrix contains empty rows")')
1363  end if
1364 
1365  ! write(20,'("m=",i6,";")')r%m
1366  ! write(20,'("n=",i6,";")')r%n
1367  ! write(20,'("nz=",i10,";")')r%nz
1368  ! call qrm_prnt_array(r%adata%rperm,'rperm',20)
1369  ! call qrm_prnt_array(r%adata%cperm,'cperm',20)
1370  ! call qrm_prnt_array(r%irn,'irn',20)
1371  ! call qrm_prnt_array(r%jcn,'jcn',20)
1372  ! call qrm_prnt_array(r%val,'val',15, 20)
1373 
1374  end subroutine _qrm_get_r
1375 
1376 end module _qrm_spmat_mod
Generif interface for the ::_qrm_spmat_destroy routine.
rcntl
Definition: dqrm_mumps.h:101
Generic interface for the qrm_adealloc_i, qrm_adealloc_2i, qrm_adealloc_s, qrm_adealloc_2s, qrm_adealloc_3s, qrm_adealloc_d, qrm_adealloc_2d, qrm_adealloc_3d, qrm_adealloc_c, qrm_adealloc_2c, qrm_adealloc_3c, qrm_adealloc_z, qrm_adealloc_2z, qrm_adealloc_3z, routines.
subroutine _qrm_csc_to_csr(in_mat, out_mat, values)
This subroutine converts a CSC matrix into a CSR matrix. Optionally the values may be ignored (this c...
Generif interface for the ::_qrm_cntl_init routine.
This module contains the interfaces of all non-typed routines.
subroutine _qrm_cntl_init(qrm_spmat)
This subroutine initializes a qrm_spmat_type instance setting default values into the control paramet...
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
subroutine _qrm_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
subroutine _qrm_psetr(qrm_spmat, string, rval)
This subroutine is meant to set the real control parameters.
subroutine _qrm_spmat_init(qrm_spmat)
This subroutine initializes a qrm_spmat_type instance setting default values into the control paramet...
subroutine _qrm_get_r(qrm_mat, r)
icntl
Definition: dqrm_mumps.h:91
This module contains the definition of the analysis data type.
subroutine _qrm_fdata_destroy(qrm_fdata)
Destroys a _qrm_fdata_type instance.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.
Generif interface for the ::_qrm_spmat_alloc routine.
Generif interface for the ::_qrm_spmat_copy routine.
This module contains all the error management routines and data.
subroutine _qrm_coo_to_csc(in_mat, out_mat, values)
This subroutine converts a COO matrix into a CSC matrix. Optionally the values may be ignored (this c...
subroutine qrm_adata_destroy(adata)
Frees an qrm_adata_type instance.
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...
The main data type for the analysis phase.
Generic interface for the qrm_aalloc_i, qrm_aalloc_2i, qrm_aalloc_s, qrm_aalloc_2s, qrm_aalloc_3s, qrm_aalloc_d, qrm_aalloc_2d, qrm_aalloc_3d, qrm_aalloc_c, qrm_aalloc_2c, qrm_aalloc_3c, qrm_aalloc_z, qrm_aalloc_2z, qrm_aalloc_3z, routines.
Definition: qrm_mem_mod.F90:78
subroutine _qrm_pseti(qrm_spmat, string, ival)
This subroutine is meant to set the integer control parameters.
integer, parameter qrm_abort_
Possible actions to be performed upon detection of an error.
Generic interface for the qrm_pdealloc_i, qrm_pdealloc_2i, qrm_pdealloc_s, qrm_pdealloc_2s, qrm_pdealloc_d, qrm_pdealloc_2d, qrm_pdealloc_c, qrm_pdealloc_2c, qrm_pdealloc_z, qrm_pdealloc_2z, routines.
Definition: qrm_mem_mod.F90:98
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...
Generic interface for the qrm_prealloc_i qrm_prealloc_s qrm_prealloc_d qrm_prealloc_c qrm_prealloc_z...
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
Generif interface for the ::_qrm_spmat_convert routine.
This type defines the data structure used to store a matrix.
gstats
Definition: dqrm_mumps.h:112
The data structure meant to store all the results of the factorization phase.
Generif interface for the ::_qrm_pseti, ::_qrm_psetr and.
Generif interface for the ::_qrm_spmat_alloc routine.
This module contains the definition of the basic sparse matrix type and of the associated methods...
subroutine _qrm_pgetr(qrm_spmat, string, rval)
Gets the values of a real control parameter. This is the dual of the ::_qrm_psetr routine; the parame...
subroutine _qrm_check_spmat(qrm_spmat, op)
Check the compatibility and correctness of icntl and rcntl parameters.
This module contains the definition of all the data related to the factorization phase.
Generif interface for the ::_qrm_pgeti, ::_qrm_pgetr and.
Generic interface for the qrm_palloc_i, qrm_palloc_2i, qrm_palloc_s, qrm_palloc_2s, qrm_palloc_d, qrm_palloc_2d, qrm_palloc_c, qrm_palloc_2c, qrm_palloc_z, qrm_palloc_2z, routines.
Definition: qrm_mem_mod.F90:57
int i
Definition: secs.c:40
subroutine _qrm_spmat_convert(in_mat, out_mat, fmt, values)
This subroutine converts an input matrix into a different storage format. Optionally the values may b...
subroutine _qrm_spmat_alloc(qrm_spmat, nz, m, n, fmt)
This subroutine allocates memory for a sparse matrix.
subroutine _qrm_spmat_copy(in_mat, out_mat, values)
This subroutine makes a copy of a matrix. Optionally the values may be ignored (this comes handy duri...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
This module contains various string handling routines.
This type defines a data structure containing all the data related to a front.
Generif interface for the ::_qrm_spmat_init routine.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.