Текст подпрограммы и версий ( Фортран )
ef02r.zip 		Тексты тестовых примеров ( Фортран )
tef02r.zip
Текст подпрограммы и версий ( Си )
ef02r_c.zip 		Тексты тестовых примеров ( Си )
tef02r_c.zip
Текст подпрограммы и версий ( Паскаль )
ef02r_p.zip 		Тексты тестовых примеров ( Паскаль )
tef02r_p.zip

Подпрограмма:  EF02R

Назначение

Приближенное решение интегрального уравнения Фредгольма I рода в классе кусочно - выпуклых функций с ограниченной первой производной.

Математическое описание

Решение интегрального уравнения Фредгольма I рода 

            b
           ∫ k(x, y)  u(y) dy  =  f(x)  ,  c ≤ x ≤ d ,
          a

при априорно известной информации об ограниченности первой производной искомого решения  u = u (y) и об участках знакопостоянства его кривизны сводится на основе метода квазирешений В.К.Иванова к задаче минимизации 

                   d      b   
     min        ∫   (  ∫  k(x, y)  u(y) dy - f(x) )2 dx
    v, u (a)    c      a

на множестве кусочно - монотонных функций  v (y) = u' (y), определенных на отрезке [a, b] 

                                 y
          u(y)  =  u(a) + ∫ v(t) dt ,   a ≤ y ≤ b ,
                               a
 и удовлетворяющих ограничениям
        v1(y) ≤ v(y) ≤ v2(y) ,   a ≤ y ≤ b ,

v1 (y) и  v2 (y) - заданные на отрезке [a, b] функции. Значение  u (a) может считаться как известным, так и подлежащим определению.

Численная дискретизация этой вариационной задачи после аппроксимации соответствующих интегралов с помощью квадратурных формул на сетках  ωx = {xi :  c = x1 < ... < xi < ... < xM = d} и ωy = {yj :  a = y1 < ... < yj < ... < yN = b} приводит к задаче квадратичного программирования 

                        M            N
  (1)      min      ∑  pi  (   ∑  dj ki j uj  -  fi  )2
            v, u1    i=1         j=1

при условии, что вектоp  v = (v1, ..., vN - 1) удовлетворяет ограничениям 

                               j -1
  (2)     uj  =  u1  +   ∑   he ve ,     j = 2, ..., N
                               e=1
        
  (3)     v1 j ≤ vj ≤ v2 j ,     j = 1, ..., N-1

  (4)     mj ( vj + 1 - vj ) ≥ 0 ,     j = 1, ..., N-2 ,

где  ki j = k (xi, yj),  fi = f (xi),  dj > 0 - коэффициенты квадратурной формулы трапеции на сетке  ωy ,  pi > 0 - весовые коэффициенты,  hj = yj + 1 - yj - шаг сетки  ωy ,  v1j = v1 (yj),  v2j = v2 (yj), параметры  mj определяются условиями кусочной монотонности функции  v = u' и принимают значения 1, 0 или -1.

Численное решение задачи (1) - (4) основано на методе проекции сопряженных градиентов, который позволяет проводить итерационный процесс минимизации в конечномерном пространстве  W2, ρ1 с нормой 

                                      N-1
    || v ||  =  ( || v || Q2  +   ∑   ρj (vj - vj-1)2 )1/2 ,
                                      j=2
                       N-1
    || v || Q  =  (   ∑   qj vj2 ) ,
                       j=1

qj > 0 и  ρj ≥ 0 - заданные числа. Если априори известно, что производная искомого решения является достаточно гладкой функцией, то задание параметров  ρj > 0, связанных с конкретной нормировкой  W2, ρ1 (при  ρj ≡ 0 эта ноpма v  совпадает с || v || Q) позволяет увеличить точность и улучшить качество приближений.

Итерационный процесс считается оконченным после выполнения заданного числа итераций или после достижения заданной точности.

1.  А.Н.Тихонов, В.Я.Арсенин, Методы решения некорректных задач. M. "Hаука", 1974.
2.  В.А.Морозов, Н.Л.Гольдман, М.К.Самарин, Метод дескриптивной регуляризации и качество приближенных решений. ИФЖ, т.33, N 6, 1977.
3.  Ф.П.Васильев, Лекции по методам решения экстремальных задач. Изд - во МГУ, 1974.
4.  Н.Л.Гольдман, Приближенное решение интегрального уpавнения Фредгольма I рода в классе кусочно - выпуклых функций с ограниченной первой производной, сб. "Численный анализ на ФОРТРАНе". Изд - во МГУ, 1978.

Использование 

    SUBROUTINE  EF02R (M, N, A, F, NG, V1, V2, MU, H, P, Q,
                                            R, E, ISP, LU, U, FU, KN, B) 
   Параметры
M - заданная размерность вектоpа правой части  f (тип: целый);
N - заданная размерность вектоpа искомого решения  u = (u1, ..., uN) (тип: целый);
A - заданный вещественный двумерный массив размерности  M * N, элементы которого  A (I, J) = ki j;
F - заданный вещественный вектоp длины  M, содержащий компоненты вектоpа правой части  f;
NG - заданный признак наличия ограничений на производную искомого решения (тип: целый);
NG = 0 - если ограничения отсутствуют;
NG = 1 - если выполнены ограничения (3);
NG = 2 - если выполнены ограничения (4);
NG = 3 - если выполнены ограничения (3) и (4);
V1 - вещественный вектоp длины  N, первые  N - 1 компонент которого содержат нижние ограничения  v1 j;
V2 - вещественный вектоp длины  N, первые  N - 1 компонент которого содержат верхние ограничения  v2 j;
MU - целый вектоp длины  N, первые  N - 2 компонент которого содержат параметры  mj;
H - вещественный вектоp длины  N, первые  N - 1 компонент которого содержат шаги  hj разностной сетки  ωy на отрезке [a, b];
P - вещественный вектоp длины  M, содержащий весовые коэффициенты  pi > 0;
Q - вещественный вектоp длины  N, первые  N - 1 компонент которого содержат весовые коэффициенты  qj > 0;
R - вещественный вектоp длины  N + 1 значений параметров  ρj таких, что:  R (1) = 0.,  R (J) = ρj ≥ 0,  j = 2, ..., N - 1,  ρ (N) = 0.,  ρ (N + 1) = 0.;
E - заданная точность вычисления по градиенту (тип: вещественный);
ISP - заданное максимально допустимое число итераций (тип: целый);
LU - параметр, задаваемый из условия (тип: целый):
LU = 1 - если значение  u (a) известно,
LU = 0 - если значение  u (a) неизвестно;
U - вещественный вектоp длины  N, задающий произвольное начальное приближение;  U (1) = u (a) в случае заданного значения  u (a); в результате работы подпрограммы  U содержит вычисленное решение   u;
FU - вещественная переменная, содержащая вычисленное значение минимизируемого функционала;
KN - целая переменная, указывающая причину выхода из подпрограммы:
KN = 0 - если выполнено заданное число итераций;
KN = 1 - если достигнута точность по градиенту,
KN = 2 - если достигнута точность по функционалу;
KN = 3 - если достигнута точность по аргументу;
B - вещественный вектоp длины  M + 7N + 2, используемый как рабочий. 
   Версии:  нет 
   Вызываемые подпрограммы
EF01R - приближенное решение интегрального уравнения Фредгольма I рода в классе ограниченных кусочно - монотонных функций. 
   Замечания по использованию
  1. 

Kpоме искомого решения  u подпрограмма позволяет получить приближение для его производной: компоненты вычисленного вектоpа  v = (v1, ..., vN - 1) расположены в  B (M + 1), ..., B (M + N - 1).

  2. 

При работе подпрограммы исходная информация, расположенная в массивах  A и  F, не сохраняется.

  3.  Справедливы замечания 3) - 5) к использованию подпрограммы EF01R.

Пример использования 

       DIMENSION  A(11, 11), F(11), V1(11), V2(11), MU(11), H(11), 
      *                        P(11), Q(11), R(12), U(11), B(90)
       DATA  V1 /10*0.0/, V2 /10*1.1/, MU /5*1, 4*-1/
       DATA  H /10*0.1/, P /11*1.0/, Q /0.05, 9*0.1, 0.05/, R /12*0.0/
       M = 11
       N = 11
       NG = 3
       LU = 1
       E = 1.0E - 7
       ISP = 20
       DO 3  I = 1, M
       DO 3  J = 1, N
       IF(I-J) 2, 1, 2
    1 A(I, J) = 1. / Q(I)
       GO TO 3
    2 A(I, J) = 0.
    3 CONTINUE
       U(1) = -11. / 24.
       DO 4  J = 2, N
    4 U(J) = 0.
       Y = -0.5
       SF = 0.025
       DO 5  I = 1, M
       S = SF + SIN(123.456*I + 789.)
       F(I) = Y - (Y**3) / 3. + S
    5 Y = Y + H(I)
       CALL  EF02R (M, N, A, F, NG, V1, V2, MU, H, P, Q, R, E,
      *                        ISP, LU, U, FU, KN, B)

Результат:

       U  =  ( -4.58E-1,   -4.11E-1,   -3.16E-1,   -2.21E-1,   -1.26E-1, 
                  -3.13E-2,    7.87E-2,    1.73E-1,    2.68E-1,    3.64E-1,   4.58E-1 )