xref: /phasta/phSolver/compressible/rstat.f (revision 595995161822a203c8467e0e4a253d7bd7d6df32)
1        subroutine rstat (res, ilwork)
2c
3c----------------------------------------------------------------------
4c
5c This subroutine calculates the statistics of the residual.
6c
7c input:
8c  res   (nshg,nflow)   : preconditioned residual
9c
10c output:
11c  The time step, cpu-time and entropy-norm of the residual
12c     are printed in the file HISTOR.DAT.
13c
14c
15c Zdenek Johan, Winter 1991.  (Fortran 90)
16c----------------------------------------------------------------------
17c
18        include "common.h"
19        include "mpif.h"
20        include "auxmpi.h"
21c
22        dimension res(nshg,nflow)
23        dimension rtmp(nshg), nrsmax(1), ilwork(nlwork)
24        dimension Forin(4), Forout(4)
25!SCATTER        dimension irecvcount(numpe), resvec(numpe)
26c        integer TMRC
27
28
29        real*8  ftots(3,0:MAXSURF),ftot(3),spmasstot(0:MAXSURF),spmasss
30
31	ttim(68) = ttim(68) - secs(0.0)
32
33        if (numpe == 1) nshgt=nshg   ! global = this processor
34c
35c incompressible style data from flx surface
36c
37      if (numpe > 1) then
38         call MPI_ALLREDUCE (flxID(2,isrfIM), spmasss,1,
39     &        MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD,ierr)
40         call MPI_ALLREDUCE (flxID(1,isrfIM), Atots,1,
41     &        MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD,ierr)
42         call MPI_ALLREDUCE (flxID(3,:), Ftots(1,:),MAXSURF+1,
43     &        MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD,ierr)
44         call MPI_ALLREDUCE (flxID(4,:), Ftots(2,:),MAXSURF+1,
45     &        MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD,ierr)
46         call MPI_ALLREDUCE (flxID(5,:), Ftots(3,:),MAXSURF+1,
47     &        MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD,ierr)
48         call MPI_ALLREDUCE (flxID(2,:), spmasstot(:),MAXSURF+1,
49     &        MPI_DOUBLE_PRECISION,MPI_SUM, MPI_COMM_WORLD,ierr)
50      else
51         Ftots=flxID(3:5,:)
52         Atots=flxID(1,isrfIM)
53         spmasss=flxID(2,isrfIM)
54         spmasstot(:)=flxID(2,:)
55      endif
56	if(myrank.eq.0) then
57      write(44,1000)lstep+1,(spmasstot(j),j=1,5)
58      call flush(44)
59	endif
60      ftot(1)=sum(Ftots(1,0:MAXSURF))
61      ftot(2)=sum(Ftots(2,0:MAXSURF))
62      ftot(3)=sum(Ftots(3,0:MAXSURF))
63c
64c end of incompressible style
65c
66c
67c.... -------------------->  Aerodynamic Forces  <----------------------
68c
69c.... output the forces and the heat flux
70c
71        if (iter .eq. nitr) then
72          Forin = (/ Force(1), Force(2), Force(3), HFlux /)
73          if (numpe > 1) then
74          call MPI_REDUCE (Forin(1), Forout(1), 4, MPI_DOUBLE_PRECISION,
75     &                                   MPI_SUM, master,
76     &                                   MPI_COMM_WORLD,ierr)
77          endif
78          Force = Forout(1:3)
79          HFlux = Forout(4)
80          if (myrank .eq. master) then
81             write (iforce,1000) lstep+1, (Force(i), i=1,nsd), HFlux,
82     &                           spmasss
83             call flush(iforce)
84          endif
85        endif
86
87c
88c.... ----------------------->  Convergence  <-------------------------
89c
90c.... compute the maximum residual and the corresponding node number
91c
92        rtmp = zero
93        do i = 1, nflow
94          rtmp = rtmp + res(:,i)**2
95        enddo
96
97        call sumgat (rtmp, 1, resnrm, ilwork)
98
99        resmaxl = maxval(rtmp)
100
101        irecvcount = 1
102        resvec = resmaxl
103        if (numpe > 1) then
104        call MPI_ALLREDUCE (resvec, resmax, irecvcount,
105     &                    MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD,
106     &                    ierr)
107c        call MPI_REDUCE_SCATTER (resvec, resmax, irecvcount,
108c     &                    MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD,
109c     &                    ierr)
110        else
111          resmax=resmaxl
112        endif
113        nrsmax = maxloc(rtmp)
114c
115c.... correct the residuals
116c
117        if (loctim(itseq) .eq. 0) then
118          resnrm = resnrm
119          resmax = resmax
120        else
121          resnrm = resnrm
122          resmax = resmax
123        endif
124c
125c.... approximate the number of entries
126c
127        totres = resnrm / float(nshgt)
128        totres = sqrt(totres)
129        resmax = sqrt(resmax)
130        if (resfrt .eq. zero) resfrt = totres
131        jtotrs = int  ( 10.d0 * log10 ( totres / resfrt ) )
132        jresmx = int  ( 10.d0 * log10 ( resmax / totres ) )
133c
134c.... get the CPU-time
135c
136        rsec=TMRC()
137        cputme = (rsec-ttim(100))
138c
139c.... output the result
140c
141        if (myrank .eq. master) then
142          print 2000,        lstep+1, cputme, totres, jtotrs, nrsmax,
143     &                     jresmx, lGMRES,  iKs, ntotGM
144          write (ihist,2000) lstep+1, cputme, totres, jtotrs, nrsmax,
145     &                     jresmx, lGMRES,  iKs, ntotGM
146          call flush(ihist)
147        endif
148	ttim(68) = ttim(68) + secs(0.0)
149
150c
151c.... return
152c
153        return
154c
1551000    format(1p,i6,5e13.5)
1562000    format(1p,i6,e10.3,e10.3,3x,'(',i4,')',3x,'<',i6,'|',i4,'>',
157     &         ' [',i3,'-',i3,']',i10)
158c
159        end
160        subroutine rstatSclr (rest, ilwork,lgmrest,ikst)
161c
162c----------------------------------------------------------------------
163c
164c This subroutine calculates the statistics of the residual.
165c
166c input:
167c  rest   (nshg)   : preconditioned residual
168c
169c output:
170c  The time step, cpu-time and entropy-norm of the residual
171c     are printed in the file HISTOR.DAT.
172c
173c
174c Zdenek Johan, Winter 1991.  (Fortran 90)
175c----------------------------------------------------------------------
176c
177        include "common.h"
178        include "mpif.h"
179        include "auxmpi.h"
180c
181        dimension rest(nshg)
182        dimension rtmp(nshg), nrsmax(1), ilwork(nlwork)
183!SCATTER        dimension irecvcount(numpe), resvec(numpe)
184c        integer TMRC
185
186	ttim(68) = ttim(68) - secs(0.0)
187        if (numpe == 1) nshgt=nshg   ! global = this processor
188c
189c.... ----------------------->  Convergence  <-------------------------
190c
191c.... compute the maximum residual and the corresponding node number
192c
193        rtmp = zero
194        rtmp = rtmp + rest**2
195
196        call sumgat (rtmp, 1, resnrm, ilwork)
197
198        resmaxl = maxval(rtmp)
199
200continue on
201
202        irecvcount = 1
203        resvec = resmaxl
204        if (numpe > 1) then
205        call MPI_ALLREDUCE (resvec, resmax, irecvcount,
206     &                    MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD,
207     &                    ierr)
208c        call MPI_REDUCE_SCATTER (resvec, resmax, irecvcount,
209c     &                    MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD,
210c     &                    ierr)
211        else
212          resmax=resmaxl
213        endif
214        nrsmax = maxloc(rtmp)
215c
216c.... correct the residuals
217c
218        if (loctim(itseq) .eq. 0) then
219          resnrm = resnrm
220          resmax = resmax
221        else
222          resnrm = resnrm
223          resmax = resmax
224        endif
225c
226c.... approximate the number of entries
227c
228        totres = resnrm / float(nshgt)
229        totres = sqrt(totres)
230        resmax = sqrt(resmax)
231        if (resfrt .eq. zero) resfrt = totres
232        jtotrs = int  ( 10.d0 * log10 ( totres / resfrt ) )
233        jresmx = int  ( 10.d0 * log10 ( resmax / totres ) )
234c
235c.... get the CPU-time
236c
237        rsec=TMRC()
238        cputme = (rsec-ttim(100))
239c
240c.... output the result
241c
242        if (myrank .eq. master) then
243          print 2000,        lstep+1, cputme, totres, jtotrs, nrsmax,
244     &                     jresmx, lgmrest,  iKst, ntotGM
245          write (ihist,2000) lstep+1, cputme, totres, jtotrs, nrsmax,
246     &                     jresmx, lgmrest,  iKst, ntotGM
247          call flush(ihist)
248        endif
249        if(totres.gt.1.0e-9) istop=istop-1
250
251	ttim(68) = ttim(68) + secs(0.0)
252
253c
254c.... return
255c
256        return
257c
2581000    format(1p,i6,4e13.5)
2592000    format(1p,i6,e10.3,e10.3,3x,'(',i4,')',3x,'<',i6,'|',i4,'>',
260     &         ' [',i3,'-',i3,']',i10)
261c
262        end
263