Учебное пособие C#. Алгоритмы и структуры данных н. А. Тюкачев, В. Г. Хлебостроев издание третье, стереотипное 1 / 23

Листинг 6.1. Метод Гаусса 
public static void Gauss(float[,] A, float[] B,
2 / 23

---

210 
ref float[] X) 
{ 
int Count = B.Length;
// 
решение системы: прямой ход 
for (int i = 0; i <= Count - 2; i++)
for (int j = i + 1; j <= Count-1; j++)
{ 
for (int k=i+1; k <= Count-1; k++)
A[k,j] += -A[i,j]*A[k,i]/A[i,i]; 
B[j] += -B[i]*A[i,j] / A[i, i]; 
} 
// 
решение системы: обратный ход 
for (int j = Count-1; j >= 0; j--)
{ 
X[j] = B[j]; 
for (int k = j + 1; k <= Count-1; k++)
X[j] += -A[k, j] * X[k]; 
X[j] /= A[j, j]; 
} 
} 
6.2. П
РИБЛИЖЕННОЕ ВЫЧИСЛЕНИЕ ПРОИЗВОДНЫХ
 
Задача численного дифференцирования состоит в приближенном вы-
числении производных функции u(x) по значениям этой функции, заданным в 
конечном числе точек. Пусть на [a, b] введена сетка 
{
}
, 0,1,
,
h
i
w
x
a ih i
N hN
b a
=
= +
=
…
= −
и определены значения u
i
=u(x
i
) функции u(x) в точках сетки. В качестве при-
ближенного значения u
′(x
i
) можно, например, взять любое из следующих 
разностных отношений: 
3 / 23

---

211 
1
1
1
1
1
,
,
,
,
,
.
2
i
i
i
i
i
i
x i
x i
x i
u
u
u
u
u
u
u
u
u
h
h
h
−
+
−
+
−
−
−
−
=
=
=
В качестве примера рассмотрим функцию 
( )
3
f x
x
= , 
производная которой равна 
( )
2
3
f x
x
′
=
. 
Код программы представлен в листинге 6.2. 
Листинг 6.2. Приближенное вычисление производной 
float a = 0; 
float b = 3; 
int m = 50; 
float dx = (b - a) / m; 
float x = 0, y = 0, yt = 0;
for (int i = 1; i <= m; i++) 
{ 
x += dx; 
y = (F3(x) - F3(x - dx)) / dx; // 
левая разность 
yt = (F3(x + dx) - F3(x - dx)) / (2 * dx);
g.DrawEllipse(Pens.Black,II(x)-2,JJ(y)-2, 4, 4); 
// 
центральная разность 
g.DrawRectangle(Pens.Black,II(x)-2,JJ(yt)-2,4,4); 
} 
Результаты расчетов представлены на рисунке 6.1. 
На рисунке 6.1 квадратами представлены приближенные значения про-
изводной, вычисленной через центральные разности. Видно, что эти значения 
совпадают с точными значениями производной, обозначенными в соответст-
вующих точках кружками. 
 
4 / 23

---

212 
Рис. 6.1. Результаты вычисления производной функции x
3
6.3. П
РИБЛИЖЕННОЕ ВЫЧИСЛЕНИЕ ИНТЕГРАЛОВ
 
Все приближенные методы вычисления определенных интегралов 
( )
b
a
I
f x dx
=

(6.18) 
основаны на замене интеграла конечной суммой 
( )
0
,
N
n
k
k
k
I
c f x
=
=

где c
k
– числовые коэффициенты и x
k
– точки отрезка [a, b], k = 0, 1, …, N. 
Приближенное равенство 
( )
( )
0
b
N
k
k
k
a
f x dx
c f x
=
≈


(6.19) 
5 / 23

---

213 
называется квадратурной формулой, а сумма вида (6.19) – квадратурной 
суммой. Точки x
k
называются узлами квадратурной формулы, а числа c
k
– 
коэффициентами квадратурной формулы. 
Введем на [a, b] равномерную сетку с шагом h, т. е. множество точек 
{
}
,
0,1 ,
,
h
i
w
x
a ih i
N hN
b a
=
= +
=
…
= −
, 
и представим интеграл (1) в виде суммы интегралов по частичным отрезкам: 
( )
( )
1
1
.
i
i
x
b
N
i
a
x
f x dx
f x dx
−
=
=



(6.20) 
Для построения формулы численного интегрирования на всем отрезке 
[a, b] достаточно построить квадратурную формулу для интеграла 
( )
1
i
i
x
x
f x dx
−

(6.21) 
на частичном отрезке [x
i-1
, x
i
] и воспользоваться свойством (6.20). 
6.3.1. Ф
ОРМУЛА ПРЯМОУГОЛЬНИКОВ
 
Заменим интеграл (6.18) выражением f(x
i-1/2
)h, где x
i-1/2
=x
i
-h/2. Геомет-
рически такая замена означает, что площадь криволинейной трапеции ABCD 
заменяется площадью прямоугольника ABCD (рис. 6.2). Тогда получаем 
формулу 
( )
1
1
2
,
i
i
x
i
x
f x dx
f x
h
−
−


≈





(6.22) 
которая называется формулой прямоугольников на частичном отрезке 
[x
i-1
, x
i
]. 
Рис. 6.2. Геометрический смысл формулы прямоугольников 
6 / 23

---

214 
Суммируя равенства (6.22) по i от 1 до N, получаем составную формулу 
прямоугольников 
( )
( )
1
1
2
.
b
N
i
i
a
f x dx
f x
h
−
=
=


(6.23) 
На рисунке 6.3 показан результат вычисления интеграла по формуле 
прямоугольников. 
Рис. 6.3. Вычисление интеграла по формуле прямоугольников 
6.3.2. Ф
ОРМУЛА ТРАПЕЦИЙ 
 
На частичном отрезке эта формула имеет вид 
( )
( )
1
1
(
.
2
i
i
x
i
i
x
f x
f x
f x dx
h
−
−
+
=

(6.24) 
Она получается путем замены подынтегральной функции f(x) интерпо-
ляционным многочленом первой степени, построенным по узлам x
i-1
, x
i
, то 
есть функцией 
( )
(
) ( ) (
) ( )
1,
1
1
1
(
).
i
i
i
i
i
L
x
x x
f x
x x f x
h
−
−
=
−
− −
Составная формула трапеций имеет вид 
( )
( ) ( )
(
)
1
0
1
1
1
0,5
0,5
,
2
b
N
i
i
N
N
i
a
f x
f x
f x dx
h h
f
f
f
f
−
−
=
+
≈
=
+ +…+
+


(6.25) 
7 / 23

---

215 
где
( )
, 0,1, , ,
i
i
f
f x
i
N hN
b a
=
=
…
= −
. 
На рисунке 6.4 показан результат вычисления интеграла по формуле 
трапеций. 
Рис. 6.4. Вычисление интеграла по формуле трапеций 
6.3.3.
 
Ф
ОРМУЛА 
С
ИМПСОНА
 
При аппроксимации интеграла (6.21) заменим функцию f(x) параболой, 
проходящей через точки (x
j
,f(x
j
)), j=i-1,i-0,5,i, т.е. представим приближенно 
f(x) в виде 
( )
( )
[
]
2,
1
, ,
,
i
i
i
f x
L
x
x
x
x
−
≈
∈
где 
( )
2,i
L
x
– интерполяционный многочлен Лагранжа второй степени, 
( )
(
)
(
)(
)
(
)
2,
1
1
1
1
1
2
2
2
2
2
2
.
i
i
i
i
i
i
i
i
i
i
L
x
x x
x x f
x x
x x f
x x
x x
f
h
−
−
−
−
−






=
−
−
−
−
−
+ −
−














(6.26) 
Проводя интегрирование, получим 
( )
(
)
1
2,
1
1
1
2
4
,
.
6
i
i
x
i
i
i
i
i
i
x
h
L
x dx
f
f
f
h x
x
−
−
−
−
=
+
+
= −

Таким образом, приходим к приближенному равенству 
8 / 23

---

216 
( )
(
)
1
2,
1
1
2
4
6
i
i
x
i
i
i
i
x
h
L
x dx
f
f
f
−
−
−
≈
+
+

, 
которое называется формулой Симпсона или формулой парабол. 
На всем отрезке [a,b] формула Симпсона имеет вид 
( )
(
)
(
)
(
)
(
)
1
1
1
2
0
1
2
1
1
3
1/ 2
2
2
4
6
2
4
.
6
b
N
i
i
i
i
a
N
N
N
h
f x dx
f
f
f
h
f
f
f
f
f
f
f
f
−
−
=
−
−
≈
+
+
=
=
+
+
+ +…+
+
+
+…+


Чтобы не использовать дробных индексов, можно обозначить 
( )
0,5 , ,
0,1, , 2 ,
i
i
i
x
a
hi f
f x
i
N hN
b a
= +
=
=
…
= −
и записать формулу Симпсона в виде 
( )
(
) (
)
(
)
0
2
2
4
2
2
1
3
2
1
2
4
.
6
b
N
N
N
a
b a
f x dx
f
f
f
f
f
f
f
f
N
−
−
−
≈
+
+
+ +…+
+
+ +…+

(6.27) 
На рисунке 6.5 показан результат вычисления интеграла по формуле 
Симпсона. 
Рис. 6.5. Вычисление интеграла по формулам Симпсона 
Все три приближенных способа представлены функциями MethodRec-
tangle(), MethodTrapezes() и MethodSimpson() в листинге 6.3. 
9 / 23

---

217 
Листинг 6.3. Приближенное вычисление интеграла 
static float F(float x) 
{ 
return 3 * x * x; 
} 
static float F1(float x) 
{ 
return x * x * x; 
} 
static float MethodRectangle(float a, float b) 
{ 
int n = 200; 
float dx = (b - a) / n; 
float result = 0; 
for (int i = 1; i <= n; i++) 
result += dx * F(a+i*dx); 
return result; 
} 
static float MethodTrapezes(float a, float b) 
{ 
int n = 200; 
float dx = (b - a) / n; 
float result = (F(a)+F(b))*dx/2; 
for (int i = 1+1; i <= n-1; i++) 
result += dx * F(a + i * dx); 
10 / 23

---

218 
return result; 
} 
static float MethodSimpson(float a, float b) 
{ 
int n = 200; 
float dx = (b - a) / n; 
int[] k = new int[] { 2, 4 }; 
float result = (F(a) + F(b)) * dx / 3; 
for (int i = 1 + 1; i <= n - 1; i++) 
result += k[i % 2]*dx*F(a + i * dx)/3; 
return result; 
} 
static void Main(string[] args) 
{ 
float a = 0; 
float b = 1; 
Console.WriteLine("
Точное={0}",F1(b)-F1(a)); 
float s = MethodRectangle(a, b); 
Console.WriteLine("
Метод прямоугольников"); 
Console.WriteLine("s = {0}", s); 
Console.WriteLine(""); 
s = MethodTrapezes(a, b); 
Console.WriteLine("
Метод трапеций"); 
Console.WriteLine("s = {0}", s); 
Console.WriteLine(""); 
11 / 23

---

219 
s = MethodSimpson(a, b); 
Console.WriteLine("
Метод Симпсона"); 
Console.WriteLine("s = {0}", s); 
Console.WriteLine(""); 
Console.ReadKey(); 
} 
В качестве примера рассмотрим интеграл 
1
2
3 1
0
0
3
|
1
x dx
x
=
=

. 
Полученные результаты показывают следующее: