ETS  \$Id: Doxyfile 2162 2020-02-26 14:16:09Z g2dpc $
 All Classes Files Functions Variables Pages
B_sol.f
Go to the documentation of this file.
1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2  subroutine solint(imov)
3 
4  include 'double.inc'
5  include 'dim.inc'
6  include 'compol.inc'
7 
8  common /comwrp/ rsp1,p1,ip1
9  real*8 zw(neqp),rsp1(nspp)
10  dimension wght(nrp)
11  integer p1(neqp),ip1(neqp),isp1(nspp),ipath,flag,esp
12  equivalence(rsp1(1),isp1(1))
13  common
14  * /c_kpr/kpr
15 
16  ! write(6,*) 'solint:enter'
17 
18  if(imov.eq.0) then
19  if(itin.eq.1) then
20  ipath=1
21  call odrvd(neqpla,ia,ja,a,p1,ip1,nspp,isp1,1,flag)
22  else
23  ipath=3
24  endif
25  call sdrvd(neqpla,p1,ip1,ia,ja,a,right,zw,nspp,
26  * isp1,rsp1,esp,ipath,flag)
27  go to 1000
28  endif
29 
30  if( itin/nitdel*nitdel+nitbeg .eq. itin
31  * .OR. itin.lt.nitbeg
32  * .OR. errm.gt.7.7d-2 ) then
33 
34 
35  call odrvd(neqpla,ia,ja,a,p1,ip1,nspp,isp1,1,flag)
36 
37  !do 10 i=1,neqpla
38  ! ip1(i)=i
39  ! p1(i)=i
40 ! 10 continue
41 
42  !write(6,*) 'odrv flag=',flag
43 
44  ipath=1
45 
46  20 continue
47 
48  call sdrvd(neqpla,p1,ip1,ia,ja,a,right,zw,nspp,
49  * isp1,rsp1,esp,ipath,flag)
50 
51 c do 860 i=1,neq
52 c write(6,*) 'zw(i) i',i,zw(i)
53 c860 continue
54 
55  if(kpr.eq.1) then
56  write(6,*) 'sdrv: flag,esp',flag,esp
57  endif
58  else
59 
60  call solbit(zw)
61 
62  endif
63 
64  1000 continue
65 
66 c raspakovka reshenia
67 
68  errpss=0.d0
69  do i=1,iplas-1
70  if(ngav.eq.0) then
71  wght(i)=0.5d0
72  !wght(i)=0.25d0
73  else
74  wght(i)=1.d0/(2.0d0+sqrt(1.d0-psia(i))*q(i)/q(1))
75  if(iswtch.eq.1) wght(i)=0.5d0
76  endif
77 
78  if(iter.eq.1) wght(i)=1.0d0
79  enddo
80 
81 
82 
83  do 500 i=1,iplas-1
84 
85  do 500 j=2,nt1
86 
87  ieq=numlin(i,j,nr,nt)
88 
89  delpsi=dabs(psi(i,j)-zw(ieq))
90  errpss=dmax1(errpss,delpsi)
91 
92  psi(i,j)=zw(ieq)*wght(i)+(1.d0-wght(i))*psi(i,j)
93 
94  500 continue
95 
96  do 600 i=1,iplas1
97 
98  psi(i,1)=psi(i,nt1)
99  psi(i,nt)=psi(i,2)
100 
101  600 continue
102 
103  !write(6,*) 'solint:errpss',errpss
104 
105  return
106  end
107 
108 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
109 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
110 
111 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
112 
113  subroutine solbit(zw)
114 
115  include 'double.inc'
116  include 'dim.inc'
117  include 'compol.inc'
118 
119  common /comwrp/ rsp1,p1,ip1
120  real*8 zw(neqp),zyy(neqp),rsp1(nspp),wpp(neqp),wzz(neqp),
121  * wrr(neqp),zuu(neqp)
122  integer p1(neqp),ip1(neqp),isp1(nspp),ipath,flag,esp
123  equivalence(rsp1(1),isp1(1))
124 
125  !write(6,*) 'solbit:enter'
126 
127  itmax=10
128 
129  relerr= 1.d-8
130 
131 
132  do i=1,nr1
133  do j=2,nt1
134 
135  ieq=numlin(i,j,nr,nt)
136 
137  zw(ieq)=psi(i,j)
138 
139  enddo
140  enddo
141 
142  do il=1,neqpla
143  wrr(il)=right(il)
144  enddo
145 
146  do il=1,neqpla
147 
148  i1=ia(il)
149  i2=ia(il+1)-1
150 
151  znes=0.d0
152 
153  do im=i1,i2
154 
155  ic=ja(im)
156 
157  znes=znes+dapp(im)*zw(ic)
158 
159  enddo
160 
161  zyy(il)=wrr(il)-znes
162 
163  enddo
164 
165  call sdrvd(neqpla,p1,ip1,ia,ja,app0,zyy,zw,nspp,
166  * isp1,rsp1,esp,3,flag)
167 
168  !if(itin/2*2.ne.itin)
169 
170  return
171 
172  itk=0
173 
174  ido=0
175 
176 
177  77 call dpcgrc(ido,neqpla,zw,wpp,wrr,wzz,relerr,itmax)
178 
179  !write(6,*) 'ido',ido
180  ! pause 'soleit '
181 
182  if(ido.eq.1) then
183 
184  do 10 il=1,neqpla
185 
186  i1=ia(il)
187  i2=ia(il+1)-1
188 
189  znes=0.d0
190 
191  do 100 im=i1,i2
192 
193  ic=ja(im)
194 
195  znes=znes+a(im)*wpp(ic)
196 
197  100 continue
198 
199  wzz(il)=znes
200 
201  10 continue
202 
203  go to 77
204 
205  elseif(ido.eq.2) then
206 
207  call sdrvd(neqpla,p1,ip1,ia,ja,app0,wrr,wzz,nspp,
208  * isp1,rsp1,esp,3,flag)
209 
210  ! write(6,*) 'sdrv: flag,esp',flag,esp
211 
212  if(itk.eq.1) go to 1056
213 
214  do il=1,neqpla
215 
216  i1=ia(il)
217  i2=ia(il+1)-1
218 
219  znes=0.d0
220 
221  do im=i1,i2
222 
223  ic=ja(im)
224 
225  znes=znes+dapp(im)*wzz(ic)
226 
227  enddo
228 
229  zyy(il)=znes
230 
231  enddo
232 
233  call sdrvd(neqpla,p1,ip1,ia,ja,app0,zyy,zuu,nspp,
234  * isp1,rsp1,esp,3,flag)
235 
236  ! write(6,*) 'sdrv: flag,esp',flag,esp
237 
238  do il=1,neqpla
239 
240  wzz(il)=wzz(il)-zuu(il)
241 
242  enddo
243 
244  1056 continue
245 
246 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
247 ! do il=1,neq
248 !
249 ! i1=ia(il)
250 ! i2=ia(il+1)-1
251 !
252 ! znes=0.d0
253 !
254 ! do im=i1,i2
255 !
256 ! ic=ja(im)
257 !
258 ! if(ic.eq.il) then
259 ! adiag=a(im)
260 ! go to 110
261 ! endif
262 !
263 ! enddo
264 !
265 ! 110 wzz(il)=wrr(il)/adiag
266 !
267 ! enddo
268 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
269 
270  ! zermax=0.d0
271  ! znvmax=0.d0
272  ! do il=1,neq
273  !
274  ! znev=dabs(wrr(il))
275  ! znvmax=dmax1(znev,znvmax)
276 
277  ! zerr=dabs(zw(il)-zw0(il))
278  ! zermax=dmax1(zerr,zermax)
279 
280  ! enddo
281 
282  itk=itk + 1
283  ! write(6,*) 'nev,err iter ',znvmax,zermax,itk
284 
285 
286  ! if(znvmax.lt.relerr) go to 537
287 
288  go to 77
289 
290  endif
291 
292  537 continue
293 
294  !write(6,*) 'nev,err iter ',znvmax,zermax,itk
295 c write(6,*) ' **iter ',itk
296 
297  ! pause ' '
298  ! stop
299  return
300  end
301 
302 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
303 
304  subroutine nev_b(zw)
305 
306  include 'double.inc'
307  include 'dim.inc'
308  include 'compol.inc'
309 
310  real*8 zw(neqp)
311 
312  do 10 il=1,neq
313 
314  i1=ia(il)
315  i2=ia(il+1)-1
316 
317  znev=0.d0
318  znmx=0.d0
319 
320  do 100 im=i1,i2
321 
322  ic=ja(im)
323 
324  znev=znev+a(im)*zw(ic)
325 
326  ! if(il.eq.1) then
327  !
328  ! write(6,*) 'im a(im)',im,a(im)
329  ! write(6,*) 'ic zw(ic)',ic,zw(ic)
330  ! write(6,*) 'znev ' ,znev
331  ! endif
332 
333  100 continue
334  znev=right(il)-znev
335  znab=dabs(znev)
336 
337 c write(6,*) '***right znev ***',il,right(il),znev
338 
339  znmx=dmax1(znmx,znab)
340 
341  10 continue
342 
343 c write(6,*) 'nev:',znmx
344 
345  return
346  end
347 
348 
349 
350 
351 
352 
353 
354 
solbit
subroutine solbit(zw)
Definition: B_sol.f:113
odrvd
subroutine odrvd
Definition: SPAR.f:9
solint
subroutine solint(imov)
Definition: B_sol.f:2
numlin
function numlin(i, j, nr, nt)
Definition: com_sub.f:1040
dpcgrc
subroutine dpcgrc(IDO, N, X, P, R, Z, RELERR, ITMAX)
Definition: DPCGRCV.f:1
sdrvd
subroutine sdrvd
Definition: SPAR.f:1411
nev_b
subroutine nev_b(zw)
Definition: B_sol.f:304