QR_MUMPS
qrm_readmat.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 
53 subroutine _qrm_readmat(matfile, qrm_mat, fakec)
54 
55  use _qrm_spmat_mod
56  use qrm_error_mod
57 
58  implicit none
59 
60  character(len=*), intent(in) :: matfile
61  type(_qrm_spmat_type), intent(inout) :: qrm_mat
62  logical, optional :: fakec
63 
64  character(len=20) :: rep, field, symm, typ, fmt
65  integer :: m, n, nnz, ival, i, myid, info, nprocs
66  logical :: values, ifakec
67  _qrm_real :: rnds(2), re, im
68  ! error management
69  integer :: err_act
70  character(len=*), parameter :: name='qrm_read_mat'
71 
72  call qrm_err_act_save(err_act)
73 
74  __qrm_prnt_msg('("Reading Matrix: ",a20)')matfile
75 
76  rep = ''
77  field = ''
78  symm = ''
79  typ = ''
80  fmt = ''
81 
82  if(present(fakec)) then
83  ifakec = fakec
84  else
85  ifakec = .false.
86  end if
87 
88  open(4,file=matfile, status='OLD', action='READ', iostat=info)
89  if(info .gt. 0) then
90  call qrm_err_push(25, name, aed=matfile)
91  goto 9999
92  end if
93 
94  read(4,*)rep,typ,fmt,field,symm
95 
96  read(4,*)rep
97  do
98  if(rep(1:1) .ne. '%') exit
99  read(4,*)rep
100  end do
101 
102  backspace(4)
103 
104  read(4,*)m,n,nnz
105 
106  values = field .ne. 'pattern'
107 
108  qrm_mat%m = m
109  qrm_mat%n = n
110  qrm_mat%nz = nnz
111 
112  call qrm_palloc( qrm_mat%irn, qrm_mat%nz )
113  call qrm_palloc( qrm_mat%jcn, qrm_mat%nz )
114  call qrm_palloc( qrm_mat%val, qrm_mat%nz )
115  __qrm_check_ret(name,'qrm_palloc',9999)
116  if(values) then
117  do i=1, nnz
118 #if defined (cprec) || defined(zprec)
119  if(field .eq. 'complex') then
120  read(4,*)qrm_mat%irn(i), qrm_mat%jcn(i), re, im
121  qrm_mat%val(i) = cmplx(re,im,kind(_qrm_one))
122  else if((field.eq.'real') .or. (field.eq.'integer')) then
123  read(4,*)qrm_mat%irn(i), qrm_mat%jcn(i), re
124  if(ifakec) then
125  qrm_mat%val(i) = cmplx(re,re,kind(_qrm_one))
126  else
127  qrm_mat%val(i) = cmplx(re,_qrm_rzero,kind(_qrm_one))
128  end if
129  end if
130 #elif defined (sprec) || defined(dprec)
131  if(field .eq. 'complex') then
132  read(4,*)qrm_mat%irn(i), qrm_mat%jcn(i), qrm_mat%val(i), im
133  else if((field.eq.'real') .or. (field.eq.'integer')) then
134  read(4,*)qrm_mat%irn(i), qrm_mat%jcn(i), qrm_mat%val(i)
135  end if
136 #endif
137  end do
138  else
139  do i=1, nnz
140  read(4,*)qrm_mat%irn(i),qrm_mat%jcn(i)
141  end do
142 #if defined (cprec) || defined(zprec)
143  qrm_mat%val = _qrm_one
144 #elif defined (sprec) || defined(dprec)
145  qrm_mat%val = _qrm_one
146 #endif
147  end if
148 
149  close(4)
150  __qrm_prnt_msg('("Matrix read.")')
151 
152  qrm_mat%fmt= 'coo'
153 
154  call qrm_err_act_restore(err_act)
155  return
156 
157 9999 continue ! error management
158  call qrm_err_act_restore(err_act)
159  if(err_act .eq. qrm_abort_) then
160  call qrm_err_check()
161  end if
162  return
163 
164 end subroutine _qrm_readmat
165 
166 
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.
This module contains all the error management routines and data.
integer, parameter qrm_abort_
Possible actions to be performed upon detection of an error.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
This type defines the data structure used to store a matrix.
subroutine _qrm_readmat(matfile, qrm_mat, fakec)
This subroutine reads a Matrix Market matrix from a file and stores it on the host processor...
Definition: qrm_readmat.F90:54
This module contains the definition of the basic sparse matrix type and of the associated methods...
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.