GNU Octave  4.0.0 A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
ssolsy.f
Go to the documentation of this file.
1  SUBROUTINE ssolsy (WM, IWM, X, TEM)
2 C***BEGIN PROLOGUE SSOLSY
3 C***SUBSIDIARY
4 C***PURPOSE ODEPACK linear system solver.
5 C***TYPE SINGLE PRECISION (SSOLSY-S, DSOLSY-D)
6 C***AUTHOR Hindmarsh, Alan C., (LLNL)
7 C***DESCRIPTION
8 C
9 C This routine manages the solution of the linear system arising from
10 C a chord iteration. It is called if MITER .ne. 0.
11 C If MITER is 1 or 2, it calls SGETRF to accomplish this.
12 C If MITER = 3 it updates the coefficient h*EL0 in the diagonal
13 C matrix, and then computes the solution.
14 C If MITER is 4 or 5, it calls SGBTRS.
15 C Communication with SSOLSY uses the following variables:
16 C WM = real work space containing the inverse diagonal matrix if
17 C MITER = 3 and the LU decomposition of the matrix otherwise.
18 C Storage of matrix elements starts at WM(3).
19 C WM also contains the following matrix-related data:
20 C WM(1) = SQRT(UROUND) (not used here),
21 C WM(2) = HL0, the previous value of h*EL0, used if MITER = 3.
22 C IWM = integer work space containing pivot information, starting at
23 C IWM(21), if MITER is 1, 2, 4, or 5. IWM also contains band
24 C parameters ML = IWM(1) and MU = IWM(2) if MITER is 4 or 5.
25 C X = the right-hand side vector on input, and the solution vector
26 C on output, of length N.
27 C TEM = vector of work space of length N, not used in this version.
28 C IERSL = output flag (in COMMON). IERSL = 0 if no trouble occurred.
29 C IERSL = 1 if a singular matrix arose with MITER = 3.
30 C This routine also uses the COMMON variables EL0, H, MITER, and N.
31 C
33 C***ROUTINES CALLED SGBTRS, SGETRS
34 C***COMMON BLOCKS SLS001
35 C***REVISION HISTORY (YYMMDD)
36 C 791129 DATE WRITTEN
37 C 890501 Modified prologue to SLATEC/LDOC format. (FNF)
38 C 890503 Minor cosmetic changes. (FNF)
39 C 930809 Renamed to allow single/double precision versions. (ACH)
40 C 010412 Reduced size of Common block /SLS001/. (ACH)
41 C 031105 Restored 'own' variables to Common block /SLS001/, to
42 C enable interrupt/restart feature. (ACH)
43 C***END PROLOGUE SSOLSY
44 C**End
45  INTEGER IWM
46  REAL WM, X, TEM
47  dimension wm(*), iwm(*), x(*), tem(*)
48  INTEGER IOWND, IOWNS,
49  1 icf, ierpj, iersl, jcur, jstart, kflag, l,
50  2 lyh, lewt, lacor, lsavf, lwm, liwm, meth, miter,
51  3 maxord, maxcor, msbp, mxncf, n, nq, nst, nfe, nje, nqu
52  REAL ROWNS,
53  1 ccmax, el0, h, hmin, hmxi, hu, rc, tn, uround
54  COMMON /sls001/ rowns(209),
55  1 ccmax, el0, h, hmin, hmxi, hu, rc, tn, uround,
56  2 iownd(6), iowns(6),
57  3 icf, ierpj, iersl, jcur, jstart, kflag, l,
58  4 lyh, lewt, lacor, lsavf, lwm, liwm, meth, miter,
59  5 maxord, maxcor, msbp, mxncf, n, nq, nst, nfe, nje, nqu
60  INTEGER I, MEBAND, ML, MU
61  REAL DI, HL0, PHL0, R
62 C
63 C***FIRST EXECUTABLE STATEMENT SSOLSY
64  iersl = 0
65  go to(100, 100, 300, 400, 400), miter
66  100 CALL sgetrs( 'N', n, 1, wm(3), n, iwm(21), x, n, inlpck)
67  RETURN
68 C
69  300 phl0 = wm(2)
70  hl0 = h*el0
71  wm(2) = hl0
72  IF (hl0 .EQ. phl0) go to 330
73  r = hl0/phl0
74  DO 320 i = 1,n
75  di = 1.0e0 - r*(1.0e0 - 1.0e0/wm(i+2))
76  IF (abs(di) .EQ. 0.0e0) go to 390
77  320 wm(i+2) = 1.0e0/di
78  330 DO 340 i = 1,n
79  340 x(i) = wm(i+2)*x(i)
80  RETURN
81  390 iersl = 1
82  RETURN
83 C
84  400 ml = iwm(1)
85  mu = iwm(2)
86  meband = 2*ml + mu + 1
87  CALL sgbtrs( 'N', n, ml, mu, 1, wm(3), meband, iwm(21), x, n,
88  * inlpck)
89  RETURN
90 C----------------------- END OF SUBROUTINE SSOLSY ----------------------
91  END
subroutine ssolsy(WM, IWM, X, TEM)
Definition: ssolsy.f:1
std::string dimension(void) const
T abs(T x)
Definition: pr-output.cc:3062