Statistics
| Branch: | Revision:

root / synthbench / euroben-dm / mod2ci / .svn / text-base / mod2ci.f.svn-base @ 0:839f52ef7657

History | View | Annotate | Download (7.8 kB)

1
      Program mod2ci
2
! ----------------------------------------------------------------------
3
! **********************************************************************
4
! *** This program is part of the Euroben Benchmark                  ***
5
! ***                                                                ***
6
! *** Copyright: EuroBen Group p/o                                   ***
7
! ***            Utrecht University, Computational Physics Dept.     ***
8
! ***            P.O. Box 80.000                                     ***
9
! ***            3508 TA Utrecht                                     ***
10
! ***            The Netherlands                                     ***
11
! ***                                                                ***
12
! *** Author of this program: Aad van der Steen                      ***
13
! *** Date                    Summer 2005                            ***
14
! **********************************************************************
15
! ----------------------------------------------------------------------
16
! --- Version 1.0 (Distributed memory, parallel, MPI).
17

    
18
! --- Purpose of program mod2ci
19
! -----------------------------
20
! --- Solver for sparse linear systems with iterative methods and
21
!     one type of preconditioner. The systems are not actually
22
!     solved but rather a predifined number of iterations is performed
23
!     to assess the speed of the solver-preconditioner combination.
24
!     Two types of systems are considered:
25
!     1. - Systems stemming from 3-D finite difference problems
26
!          resulting in 7-banded matrices. We consider only the
27
!          symmetric type in ANOTHER program: 'mod2cr'.
28
!     2. - Finite element type, irregularly filled matrices.
29
!          They are stored in CRS format. In this program 'mod2ci' we
30
!          only look at non-symmetric matrices and their solution
31
!          methods.
32
!
33
! --- In this Program 'mod2ci' we address the irregular, FEM type
34
!     of problem. For the regular 7-band matrix-based type we have
35
!     program 'mod2cr'.
36
!
37
!     Solvers used:
38
!     -  For the symmetric banded systems in 'mod2cr': CG with ILU(0)
39
!        and polynomial preconditioning.
40
!     -  For irregular systems in 'mod2ci': RGMRES and TFQMR with
41
!        polynomial preconditioning.
42
! ----------------------------------------------------------------------
43
      Use                       numerics
44
      Use                       floptime
45
      Use                       mpi_module
46
      Implicit                  None
47

    
48
      Integer, Allocatable   :: indx(:), rowp(:)
49
      Real(l_), Allocatable  :: matvals(:), b(:), q(:), x(:)
50
      Real(l_)               :: corr, maxdf, mindf, res
51
      Real(l_)               :: frac, gtime, mflops, resnrm, time
52
      Integer, Parameter     :: maxit = 50
53
      Integer, Parameter     :: mp = 1
54
      Real(l_), Parameter    :: tol = 1.0e-10_l_
55
      Real(l_), Parameter    :: micro = 1.0e-6_l_, two = 2.0_l_,
56
     &                          perc = 1.0e2_l_
57
      Real(l_)               :: gamma(mp+1)
58
      Real(l_)               :: start_time, end_time
59
      Integer                :: i, lelmts, lrows, ncols, nelmts, nrows,
60
     &                          nrep
61
      Integer(8)             :: gflops
62
      Integer                :: rg
63
      Logical                :: ok
64
      External               :: lpolyn
65
! ----------------------------------------------------------------------
66
      Call mpistart
67
      start_time = MPI_Wtime()
68
      If ( me == 0 ) Call state( 'mod2ci  ' )
69
      Open( 1, File = 'mod2ci.in' )
70
      If ( me == 0 ) Print 1000, nodes
71
! ----------------------------------------------------------------------
72
! --- Do case with GMRES and polynomial preconditioning.
73
 
74
      If ( me == 0 ) Print 1005
75
   10 Read( 1, *, End = 20 ) ncols, nrows, nelmts, nrep
76
         Call evdist( nrows )
77
         Call bsaddr
78
         lrows  = part( nrows )
79
         lelmts = part( nelmts )     
80
         Allocate( indx(lelmts), rowp(lrows+1), matvals(lelmts),
81
     &             b(ncols), q(nrows), x(nrows) )
82
         Call getmatvec( ncols, lrows, lelmts, indx, rowp, matvals, b )
83
         Call pcoefs( mp, gamma )
84
         q = 1.0_l_
85
         flops = 0
86
         time = MPI_Wtime()
87
         Do i = 1, nrep
88
            x = 0.0_l_
89
            resnrm = -1.0_l_
90
            Call rgmres( ncols, lrows, lelmts, mp, indx, rowp, matvals,
91
     &                   q, x, b, gamma, maxit, tol, resnrm, lpolyn )
92
         End Do
93
         time = MPI_Wtime() - time
94
         Call MPI_Allreduce(flops, gflops, 1, ityp, MPI_Sum, comm, ierr)
95
         Call MPI_Allreduce( time, gtime, 1, rtyp, MPI_Max, comm, ierr )
96
         mflops = micro*Real( gflops, l_ )/time
97
         gtime   = gtime/Real( nrep, l_ )
98
         frac   = perc*Real( nelmts, l_ )/( Real( ncols, l_)*
99
     &                 Real( nrows, l_ ) )
100
         If ( me == 0 ) Print 1010, ncols, nrows, frac, gtime, mflops,
101
     &                              resnrm
102
         Deallocate( indx, rowp, matvals, b, q, x )
103
      Go To 10
104
   20 If ( me == 0 ) Print 1020
105
      If ( me == 0 ) Print 1025
106
! ----------------------------------------------------------------------
107
! --- Do case with TFQMR and polynomial preconditioning.
108
 
109
      Rewind 1
110
      If ( me == 0 ) Print 1030, nodes
111
      If ( me == 0 ) Print 1005
112
   30 Read( 1, *, End = 40 ) ncols, nrows, nelmts, nrep
113
         Call evdist( nrows )
114
         Call bsaddr
115
         lrows  = part( nrows )
116
         lelmts = part( nelmts )     
117
         Allocate( indx(lelmts), rowp(lrows+1), matvals(lelmts),
118
     &             b(ncols), q(nrows), x(nrows) )
119
         Call getmatvec( ncols, lrows, lelmts, indx, rowp, matvals, b )
120
         Call pcoefs( mp, gamma )
121
         q = 1.0_l_
122
         flops = 0
123
         time = MPI_Wtime()
124
         Do i = 1, nrep
125
            x = 0.0_l_
126
            resnrm = -1.0_l_
127
            Call tfqmr( ncols, lrows, lelmts, mp, indx, rowp, matvals,
128
     &                  q, x, b, gamma, maxit, tol, resnrm, lpolyn )
129
         End Do
130
         time = MPI_Wtime() - time
131
         Call MPI_Allreduce(flops, gflops, 1, ityp, MPI_Sum, comm, ierr)
132
         Call MPI_Allreduce( time, gtime, 1, rtyp, MPI_Max, comm, ierr )
133
         mflops = micro*Real( gflops, l_ )/time
134
         gtime  = gtime/Real( nrep, l_ )
135
         frac   = perc*Real( nelmts, l_ )/( Real( ncols, l_)*
136
     &                 Real( nrows, l_ ) )
137
         If ( me == 0 ) Print 1010, ncols, nrows, frac, gtime, mflops,
138
     &                              resnrm
139
         Deallocate( indx, rowp, matvals, b, q, x )
140
      Go To 30
141
   40 If ( me == 0 ) Print 1020
142
      Call mpibye
143
      
144
      end_time = MPI_Wtime() - start_time
145
      If (me == 0) Then
146
         Write(6,22) 'Walltime: ', end_time, " s"
147
 22      Format(A,F9.3,A)
148
         Write(6,*) 'Program terminated normally'
149
      End If
150
! ----------------------------------------------------------------------
151
 1000 Format( 'Program mod2ci: Sparse iterative solver test'/
152
     &        'Non-symmetric, CRS: RGMRES with polyn. preconditioner.'/
153
     &        'No. of proc.s = ', i5 )
154
 1005 Format( '-------------------------------------------------------',
155
     &        '-----------'/
156
     &        ' #Rows | #Cols | %Fill |   Time(s)   |   Mflop/s   |   R'
157
     &        'esidue   |'/
158
     &        '-------------------------------------------------------',
159
     &        '----------|' )
160
 1010 Format( i7, '|', i7, '|', f6.2, ' |', g13.5, '|', g13.5, '|',
161
     &        g13.5, '|' )
162
 1020 Format( '-------------------------------------------------------', 
163
     &        '-----------' )
164
 1025 Format ( / )
165
 1030 Format( '-------------------------------------------------------',
166
     &        '-----------'/
167
     &        'Program mod2ci: Sparse iterative solver test'/
168
     &        'Non-symmetric, CRS: TFQMR with polyn. preconditioner.'/
169
     &        'No. of proc.s = ', i5 )
170
! ----------------------------------------------------------------------
171
      End Program mod2ci