Eötvös Quantum Utilities  v4.8.128
Providing the Horsepowers in the Quantum Realm
test_dgetPartialInv.F90
Go to the documentation of this file.
1 !#include "fintrf.h"
2 !======================================================================
3 ! Program to test the calculation of partial inverse of a real valued matrix.
4 ! Copyright (C) 2018 Peter Rakyta, Ph.D.
5 !
6 ! This program is free software: you can redistribute it and/or modify
7 ! it under the terms of the GNU General Publi! License as published by
8 ! the Free Software Foundation, either version 3 of the License, or
9 ! (at your option) any later version.
10 !
11 ! This program is distributed in the hope that it will be useful,
12 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
13 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 ! GNU General Public License for more details.
15 !
16 ! You should have received a copy of the GNU General Public License
17 ! along with this program. If not, see http://www.gnu.org/licenses/.
18 !
19 !======================================================================
20  PROGRAM matrixio
21 
22  use getpartialinv
23 
24 
25 ! Declarations
26  implicit none
27 
28 
29  integer::narg,cptArg !#of arg & counter of arg
30  character(len=100):: name !Arg name
31  logical LookForFilename, LookForinvDim
32  logical LookForOutput
33  logical fileExist
34  character(len=100) :: filename, output
35 
36 ! Array information:
37  integer*4 nonzerosA !number of nonzeros elements in A
38  integer*4 sizeAm, sizeAn !size of the matrix A (rows, and cols)
39  integer*4 sizeInv !the size of the partial inverse to be calculated
40  character(len=20) sizeInv_string !the size of the partial inverse -- string from the input
41  character(len=8) out_status
42 
43 ! Arguments for computational routine:
44  integer*4, allocatable :: ia(:), ja(:)
45 #ifdef COMP
46  COMPLEX*16, ALLOCATABLE :: a(:)
47 #else
48  REAL*8, ALLOCATABLE :: a(:)
49 #endif
50 
51 #ifdef COMP
52  COMPLEX*16, allocatable :: invA(:,:)
53 #else
54  real*8, allocatable :: invA(:,:)
55 #endif
56 !
57 
58  integer i,j
59 
60 #ifdef DEBUG
61 ! Debug parameters
62  integer*4 tmp
63  character*800 line
64 #endif
65 
66 ! https://genomeek.wordpress.com/2012/02/09/fortran-command-line-arguments/
67  !Check if any arguments are found
68  lookforfilename = .false.
69  lookforinvdim = .false.
70  lookforoutput = .false.
71  filename = ''
72  output = ''
73  sizeinv = 0
74 
75  narg=command_argument_count()
76  !Loop over the arguments
77  if(narg>0)then
78  !loop across options
79  do cptarg=1,narg
80  call get_command_argument(cptarg,name)
81  select case(adjustl(name))
82  case("--help","-h")
83  write(0,*)"This is program calculates a " , &
84  "partial inverse of a real valued sparse matirix. Usage: ./MatrixIO ", &
85  "--filename FILENAME --invDim INTEGER --Output FILENAME"
86  stop
87  case("--filename","-F", "--Filename")
88  lookforfilename = .true.
89  case("--output","-o", "--Output")
90  lookforoutput = .true.
91  case("--invDim","-I")
92  lookforinvdim = .true.
93  case default
94  if (lookforfilename) then
95  filename=adjustl(name) !assign a value to pedfile
96  lookforfilename = .false.
97  elseif (lookforoutput) then
98  output=adjustl(name) !assign a value to pedfile
99  lookforoutput = .false.
100  elseif (lookforinvdim) then
101  sizeinv_string = trim(adjustl(name))
102  READ( sizeinv_string, * ) sizeinv
103  lookforinvdim = .false.
104  else
105  write(0,*)"Option ",adjustl(name),"unknown"
106  end if
107  end select
108  end do
109  end if
110 
111  if (output.eq.'') then
112  output = filename
113  end if
114 
115 #ifdef DEBUG
116  write(0,*) 'Input file name:', trim(filename)
117  write(0,*) 'Output file name:', trim(output)
118 #endif
119 
120 ! Testing the existance of the input file containing the sparse matrix A
121  inquire(file=trim(filename),exist=fileexist)!check if it exist
122  if(.not.fileexist)then
123  write(0,*)"file ",trim(filename)," not found!"
124  stop
125  endif
126 
127 ! Opening the file to read in data
128  OPEN(unit=6, file=trim(filename),status='OLD',form='UNFORMATTED')
129 
130  READ(6,end=999,err=1000) sizean
131  READ(6,end=999,err=1000) nonzerosa
132 
133  allocate(ia(sizean+1))
134  allocate(ja(nonzerosa))
135  allocate(a(nonzerosa))
136 
137  READ(6,end=999,err=1000) ia(:)
138  READ(6,end=999,err=1000) ja(:)
139  READ(6,end=999,err=1000) a(:)
140 
141 
142 
143 ! Closing the input file
144  CLOSE(unit=6)
145 
146 #ifdef DEBUG
147  write(0,*) 'Size of the partial inverse:', sizeinv
148  write(0,*) 'The size of the input MATRIX:', sizean
149  write(0,*) 'Nonzeros in the input matrix:', nonzerosa
150  write(0,*) 'The input MATRIX:'
151  i = 0
152  do j = 1, size(a)
153  if (j >= ia(i+1)) then
154  i=i+1
155  end if
156 
157  if (j < ia(i+1)) then
158  !write(0,*) 'row', i, 'col', ja(j), 'element of A = ', a(j)
159  end if
160  end do
161 
162 #endif
163 
164 
165 
166 
167 ! Calculating the partial inverse
168  ALLOCATE(inva(sizeinv, sizeinv))
169 
170  CALL dgetpartialinv(sizean, nonzerosa, a, ia, ja, sizeinv, inva)
171 
172 ! deallocate the fortran matrix A
173  deallocate(ia)
174  deallocate(ja)
175  deallocate(a)
176 #ifdef DEBUG
177  do i=1, sizeinv
178  do j=1, sizeinv
179  write(0,*) 'row', i, 'col', j, ' element of A = ', inva(i,j)
180  end do
181  end do
182 #endif
183 
184 
185 ! Deallocating memory for the partial inverse
186 
187  if ( allocated(inva) ) then
188  DEALLOCATE(inva)
189  endif
190 
191 
192  stop
193 
194  999 write(*,'(/"End-of-file when i = ",I5)') i
195  stop
196 
197  1000 write(*,'(/"ERROR reading when i = ",I5)') i
198  stop
199  end
200 
201 
Module to calculate the partial inverse of a real/complex sparse matrix via the PARDISO libraries...
program matrixio
Definition: MatrixIO.F90:22
subroutine dgetpartialinv(sizeA, nonzerosA, a, ia, ja, sizeInv, invA, error)
Get the partial inverse of double valued real matrices.