1 SUBROUTINE sprepj (NEQ, Y, YH, NYH, EWT, FTEM, SAVF, WM, IWM,
58 REAL Y, YH, EWT, FTEM, SAVF, WM
59 dimension neq(*), y(*), yh(nyh,*), ewt(*), ftem(*), savf(*),
61 INTEGER INIT, MXSTEP, MXHNIL, NHNIL, NSLAST, CNYH,
62 1 ialth, ipup, lmax, meo, nqnyh, nslp,
63 1 icf, ierpj, iersl, jcur, jstart, kflag, l,
64 2 lyh, lewt, lacor, lsavf, lwm, liwm, meth, miter,
65 3 maxord, maxcor, msbp, mxncf, n, nq, nst, nfe, nje, nqu
66 REAL CONIT, CRATE, EL, ELCO, HOLD, RMAX, TESCO,
67 1 ccmax, el0, h, hmin, hmxi, hu, rc, tn, uround
68 COMMON /sls001/ conit, crate, el(13), elco(13,12),
69 1 hold, rmax, tesco(3,12),
70 1 ccmax, el0, h, hmin, hmxi, hu, rc, tn, uround,
71 2 init, mxstep, mxhnil, nhnil, nslast, cnyh,
72 3 ialth, ipup, lmax, meo, nqnyh, nslp,
73 3 icf, ierpj, iersl, jcur, jstart, kflag, l,
74 4 lyh, lewt, lacor, lsavf, lwm, liwm, meth, miter,
75 5 maxord, maxcor, msbp, mxncf, n, nq, nst, nfe, nje, nqu
76 INTEGER I, I1, I2, IER, II, J, J1, JJ, LENP,
77 1 mba, mband, meb1, meband, ml, ml3, mu, np1
78 REAL CON, DI, FAC, HL0, R, R0, SRUR, YI, YJ, YJJ,
86 go
to(100, 200, 300, 400, 500), miter
91 CALL jac(neq, tn, y, 0, 0, wm(3), n)
94 120 wm(i+2) = wm(i+2)*con
97 200 fac = svnorm(n, savf, ewt)
98 r0 = 1000.0e0*
abs(h)*uround*n*fac
99 IF (r0 .EQ. 0.0e0) r0 = 1.0e0
104 r =
max(srur*
abs(yj),r0/ewt(j))
107 CALL
f(neq, tn, y, ftem)
109 220 wm(i+j1) = (ftem(i) - savf(i))*fac
118 wm(j) = wm(j) + 1.0e0
121 CALL sgetrf(n, n, wm(3), n, iwm(21), ier)
122 IF (ier .NE. 0) ierpj = 1
128 310 y(i) = y(i) + r*(h*savf(i) - yh(i,2))
129 CALL
f(neq, tn, y, wm(3))
132 r0 = h*savf(i) - yh(i,2)
133 di = 0.1e0*r0 - h*(wm(i+2) - savf(i))
135 IF (
abs(r0) .LT. uround/ewt(i)) go
to 320
136 IF (
abs(di) .EQ. 0.0e0) go
to 330
137 wm(i+2) = 0.1e0*r0/di
151 CALL jac(neq, tn, y, ml, mu, wm(ml3), meband)
154 420 wm(i+2) = wm(i+2)*con
164 fac = svnorm(n, savf, ewt)
165 r0 = 1000.0e0*
abs(h)*uround*n*fac
166 IF (r0 .EQ. 0.0e0) r0 = 1.0e0
170 r =
max(srur*
abs(yi),r0/ewt(i))
172 CALL
f(neq, tn, y, ftem)
173 DO 550 jj = j,n,mband
176 r =
max(srur*
abs(yjj),r0/ewt(jj))
180 ii = jj*meb1 - ml + 2
182 540 wm(ii+i) = (ftem(i) - savf(i))*fac
189 wm(ii) = wm(ii) + 1.0e0
192 CALL sgbtrf( n, n, ml, mu, wm(3), meband, iwm(21), ier)
193 IF (ier .NE. 0) ierpj = 1
F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T const F77_REAL const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T const F77_DBLE F77_DBLE &F77_RET_T const F77_REAL F77_REAL &F77_RET_T F77_REAL F77_REAL &F77_RET_T F77_DBLE F77_DBLE &F77_RET_T const F77_DBLE const F77_DBLE * f
may be zero for pure relative error test tem the relative tolerance must be greater than or equal to
OCTAVE_EXPORT octave_value_list etc The functions then dimension(columns)
charNDArray max(char d, const charNDArray &m)
OCTAVE_EXPORT octave_value_list return the value of the option it must match the dimension of the state and the relative tolerance must also be a vector of the same length tem it must match the dimension of the state and the absolute tolerance must also be a vector of the same length The local error test applied at each integration step is xample roup abs(local error in x(i))<
subroutine sprepj(NEQ, Y, YH, NYH, EWT, FTEM, SAVF, WM, IWM, F, JAC)
charNDArray min(char d, const charNDArray &m)