18.2 Volterra Equations 795 in 818.0),but also initial-value problems for ordinary differential equations.In fact. many algorithms for ODEs have counterparts for Volterra equations. The simplest way to proceed is to solve the equation on a mesh with uniform spacing: t=a+h,i=0,1,,N, b-a h三 (18.2.2) To do so,we must choose a quadrature rule.For a uniform mesh,the simplest scheme is the trapezoidal rule,equation (4.1.11): -1 ["K(t)f()ds= K00+∑Kj+Kf月 (18.2.3) j=1 不 Thus the trapezoidal method for equation (18.2.1)is: fo go 2>、 i-1 9 (1-hKi)fi=h Kon+∑K (18.2.4) +9i i=1,.,N j=1 (For a Volterra equation of the first kind,the leading 1 on the left would be absent. and g would have opposite sign,with corresponding straightforward changes in the rest of the discussion.) Equation(18.2.4)is an explicit prescription that gives the solution in O(N2) operations.Unlike Fredholm equations,it is not necessary to solve a system of linear equations.Volterra equations thus usually involve less work than the corresponding Fredholm equations which,as we have seen,do involve the inversion of,sometimes large,linear systems. The efficiency of solving Volterra equations is somewhat counterbalanced by the fact that systems of these equations occur more frequently in practice.If we interpret equation(18.2.1)as a vector equation for the vector of m functions f(t), then the kernel K(t,s)is an m x m matrix.Equation (18.2.4)must now also be Numerica 10.621 understood as a vector equation.For each i,we have to solve the m x m set of 431 linear algebraic equations by Gaussian elimination. The routine voltra below implements this algorithm.You must supply an (outside Recipes external function that returns the kth function of the vector g(t)at the point t,and another that returns the (l)element of the matrix K(t,s)at(t,s).The routine North voltra then returns the vector f(t)at the regularly spaced points ti. #include "nrutil.h" void voltra(int n,int m,float to,float h,float *t,float **f, float (*g)(int,float),float (*ak)(int,int,float,float)) Solves a set ofm linear Volterra equations of the second kind using the extended trapezoidal rule. On input,to is the starting point of the integration and n-1 is the number of steps of size h to be taken.g(k,t)is a user-supplied external function that returns g(t),while ak(k,1,t,s) is another user-supplied external function that returns the (k,element of the matrix K(t,s). The solution is returned in f [1..m][1..n],with the corresponding abscissas in t[1..n].18.2 Volterra Equations 795 Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). in §18.0), but also initial-value problems for ordinary differential equations. In fact, many algorithms for ODEs have counterparts for Volterra equations. The simplest way to proceed is to solve the equation on a mesh with uniform spacing: ti = a + ih, i = 0, 1,...,N, h ≡ b − a N (18.2.2) To do so, we must choose a quadrature rule. For a uniform mesh, the simplest scheme is the trapezoidal rule, equation (4.1.11):  ti a K(ti, s)f(s) ds = h   1 2Ki0f0 + i−1 j=1 Kijfj + 1 2Kiifi   (18.2.3) Thus the trapezoidal method for equation (18.2.1) is: f0 = g0 (1 − 1 2hKii)fi = h   1 2Ki0f0 + i−1 j=1 Kijfj   + gi, i = 1,...,N (18.2.4) (For a Volterra equation of the first kind, the leading 1 on the left would be absent, and g would have opposite sign, with corresponding straightforward changes in the rest of the discussion.) Equation (18.2.4) is an explicit prescription that gives the solution in O(N 2) operations. Unlike Fredholm equations, it is not necessary to solve a system of linear equations. Volterra equations thus usually involve less work than the corresponding Fredholm equations which, as we have seen, do involve the inversion of, sometimes large, linear systems. The efficiency of solving Volterra equations is somewhat counterbalanced by the fact that systems of these equations occur more frequently in practice. If we interpret equation (18.2.1) as a vector equation for the vector of m functions f(t), then the kernel K(t, s) is an m × m matrix. Equation (18.2.4) must now also be understood as a vector equation. For each i, we have to solve the m × m set of linear algebraic equations by Gaussian elimination. The routine voltra below implements this algorithm. You must supply an external function that returns the kth function of the vector g(t) at the point t, and another that returns the (k,l) element of the matrix K(t, s) at (t, s). The routine voltra then returns the vector f(t) at the regularly spaced points t i. #include "nrutil.h" void voltra(int n, int m, float t0, float h, float *t, float **f, float (*g)(int, float), float (*ak)(int, int, float, float)) Solves a set of m linear Volterra equations of the second kind using the extended trapezoidal rule. On input, t0 is the starting point of the integration and n-1 is the number of steps of size h to be taken. g(k,t) is a user-supplied external function that returns gk(t), while ak(k,l,t,s) is another user-supplied external function that returns the (k, l) element of the matrix K(t, s). The solution is returned in f[1..m][1..n], with the corresponding abscissas in t[1..n]. {