Методы решения обычных дифференциальных уравнений (Лабораторная работа № 7), страница 10

C   ВЫХОДНЫЕ ФАЙЛЫ  AAAA.REZ I AAAA.GRF?

4 IF(X-H*REAL(NREZ*INTWR).LT.0.)GO TO 6

C        УВЕЛИЧЕНИЕ СЧЕТЧИКА ЗАПИСЕЙ НА ЕДИНИЦУ

NREZ=NREZ+1

C ЗАПИСЬ ИНФОРМАЦИИ ДЛЯ ИЗБРАННОГО УЗЛА В ВЫХОДНОЙ

С  ФАЙЛ AAAA.REZ

WRITE(3,103)X,(Y(I),I=1,N)

C ЗАПИСЬ ИНФОРМАЦИИ ДЛЯ ИЗБРАННОГО УЗЛА В ВЫХОДНОЙ

С  ФАЙЛ AAAA.GRF

WRITE(4,103)X,(Y(I),I=1,N)

103 FORMAT(10F8.4)

C ВЫВОД ИНФОРМ. НА ЭКРАН ДИСПЛЕЯ ДЛЯ ТЕКУЩЕГО

С  КОНТРОЛЯ РЕШЕНИЯ

WRITE(5,106)X,NREZ

106 FORMAT(1X,'X=',F8.4,'    NREZ=',I4)

C ОБРАЩЕНИЕ К ПОДПРОГРАММЕ РЕШЕНИЯ

С  ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ

6 CALL RNGKT4(N,X,H,Y,Y0,PRAV)

C         НАДО ЛИ ПРОДОЛЖАТЬ РАСЧЕТ?

9 IF(NREZ.LT.501)GO TO 10

WRITE(5,*)'DIFRK4: РАСЧЕТ ПРИОСТАНОВЛЕН ПО ПРИЧИНЕ',

*' ЗАПОЛНЕНИЯ ВЫХОДНЫХ ФАЙЛОВ (NREZ>501)'

GO TO 11

10 IF ((X-XKAN).GT.0.)GO TO 11

C  УВЕЛИЧЕНИЕ НЕЗАВИСИМОЙ ПЕРЕМЕННОЙ НА ВЕЛИЧИНУ

С  ШАГА H

X=X+H

GO TO 4 

C        ЗАКРЫТИЕ  ВЫХОДНЫХ  ФАЙЛОВ

11 CLOSE(UNIT=3)

CLOSE(UNIT=4)

WRITE(5,101) Z2,Z3

101 FORMAT(1X,'DIFRK4: РАБОТА ЗАКОНЧЕНА РЕЗУЛ. ЗАПИСАНЫ',

*' В ФАЙЛЫ ',A8,' I ',A8)                                                          422

STOP

END

SUBROUTINE RNGKT4(N,X,H,Y,Y0,FNC)

C ПОДПРОГРАММА РЕШЕНИЯ СИСТЕМЫ ДИФФЕРЕН.УРАВНЕНИЙ        

С dY/dx=F(x,Y) С НАЧАЛЬН.УСЛОВИЯМИ (ЗАДАЧА КАШИ)

C МЕТОДОМ РУНГЕ-КУТТА ЧЕТВЕРТ. ПОРЯДКА НА ОДНОМ

С ВРЕМЕННОМ ШАГЕ С ЗАДАНЫМ ШАГАМ ПО ВРЕМЕНИ

C ОБОЗНАЧЕНИЯ ПЕРЕМЕННЫХ:

C    N - КОЛИЧЕСТВО УРАВНЕНИЙ В СИСТЕМЕ (ВХОДНАЯ

С  ВЕЛИЧИНА)

C    X - ТЕКУЩЕЕ ЗНАЧЕНИЕ НЕЗАВИСИМОЙ  ПЕРЕМ. x (РАБОЧАЯ И

C    ВЫХОДНАЯ   ВЕЛИЧИНА)

C    H - ШАГ ИЗМЕНЕНИЯ НЕЗАВИСИМОЙ  ПЕРЕМЕННОЙ x

С  (ВХОДНАЯ ВЕЛИЧИНА)

C    Y0(N) - МАССИВ ЗНАЧ. ИНТЕГРИРУЕМЫХ ПЕРЕМ. Y В НАЧАЛЕ

C    КАЖДОГО  ШАГА  ПРИ ПЕРВОМ ОБРАЩЕНИИ К ПОДПР. В

C     ЭТОМ МАССИВЕ ДОЛЖНЫ НАХОДИТСЯ НАЧАЛЬНЫЕ

C   ЗНАЧЕНИЯ ИНТЕГРИРУЕМЫХ ПЕРЕМЕННЫХ  (ВХОДНЫЕ И

С  РАБОЧИЕ ВЕЛИЧИНЫ)

C    Y(N) - МАССИВ ТЕКУЩИХ ЗНАЧ. ИНТЕГРИРУЕМЫХ  ПЕРЕМ. Y В           

С  ТЕЧЕНИЕ КАЖДОГО ШАГА (РАБОЧИЙ МАСIУ)

C    YF(N) - МАССИВ ЗНАЧЕНИЙ  ИНТЕГРИРУЕМЫХ ПЕРЕМЕННЫХ

C    Y, ИСПОЛЬЗ.   В КАЧЕСТВЕ АРГУМЕНТОВ ПРИ ВЫЧИСЛЕНИИ

C    ФУНКЦИЙ F(x,Y) В ПРАВЫХ  ЧАСТЯХ УРАВНЕНИЙ  (РАБОЧИЙ

С  МАССИВ)

C    F(N) - МАССИВ ЗНАЧ. ФУНКЦИЙ F(x,Y) В ПРАВЫХ ЧАСТЯХ

C  УРАВНЕНИЙ  (РАБОЧИЙ МАССИВ).ЗНАЧЕНИЯ ФУНКЦИЙ

C   ВЫСЧИТЫВАЮТСЯ ВО ВНЕШНЕЙ  ПОДПРОГРАММЕ FNC,

С  КОТОРАЯ ЦЕЛИКОМ ПРЕДСТАВЛЕНА ВИДОМ ПРАВОЙ ЧАСТИ

C  СИСТЕМЫ УРАВНЕНИЙ I ДЛЯ КАЖДОЙ НОВОЙ СИСТЕМЫ

C    СКЛАДЫВАЕТСЯ   ПОЛЬЗОВАТЕЛЕМ ЗАНОВО

C    FNC - ИМЯ ВНЕШНЕЙ ПОДПРОГРАММЫ, В КОТОРОЙ

C  ВЫСЧИТЫВАЮТСЯ ЗНАЧЕНИЯ ФУНКЦИЙ F(x,Y) В ПРАВОЙ

C   ЧАСТИ СИСТЕМЫ  ИМЯ ЭТОЙ ПОДПР. В СПИСКЕ ФОРМАЛЬНЫХ

С  ПАРАМЕТРОВ ОПЕРАТОРА SUBROUTINE ЭТОЙ

C    ПОДПРОГРАММЫ И В СПИСКЕ ФАКТИЧЕСКИХ ПАРАМЕТРОВ

C    ОПЕРАТОРА CALL  ПРОГРАММЫ, КОТОРАЯ ЕЕ ВЫЗЫВАЕТ,

C   МОЖЕТ БЫТЬ РАЗНЫМ, В ТАКОМ СЛУЧАЕ В ПРОГРАММЕ,

С  КОТОРАЯ ЕЕ ВЫЗЫВАЕТ, ЭТО ИМЯ ДОЛЖНО

C   БЫТЬ ОПИСАНО В ОПЕРАТОРЕ EXTERNAL

C    P(4), Q(4) - РАБОЧИЕ МАССИВЫ КОЭФФИЦИЕНТОВ, ИСПОЛЬЗ.

C  ПРИ  РЕАЛИЗАЦИИ ФОРМУЛ МЕТОДОМ РУНГЕ-КУТТА

С  ЧЕТВЕРТОГО ПОРЯДКА

DIMENSION Y(*),Y0(*),YF(10),F(10),P(4),Q(4)

DATA P/0.,0.5,0.5,1./,Q/1.,2.,2.,1./  

C НАЧАЛО ЦИКЛА ПО ЧЕТЫРЕМ ТАКТАМ  РУНГЕ-КУТТА

4 DO 7 J=1,4

C ВЫЧИСЛЕНИЕ НЕЗАВИС. ПЕРЕМ. x И ИНТЕГРИРУЕМЫХ

C  ПЕРЕМЕННЫХ Y КАК  АРГУМЕНТОВ ФУНКЦИЙ F(x,Y) В ПРАВЫХ

С  ЧАСТЯХ УРАВНЕНИЙ

XF=X+P(J)*H 

DO 5 I=1,N

5 YF(I)=Y0(I)+H*F(I)*P(J)          

C ВЫЧИСЛЕНИЕ ЗНАЧ. ФУНКЦИЙ F(x,Y) В ПРАВЫХ ЧАСТЯХ

С  УРАВНЕНИЙ

CALL FNC(XF,YF,F)

C ВЫЧИСЛЕНИЕ ЗНАЧЕНИЙ ИНТЕГРИРУЕМЫХ  ПЕРЕМЕННЫХ Y

С  НА НОВОМ ТАКТЕ

DO 6 I=1,N

6 Y(I)=Y(I)+H*F(I)*Q(J)/6.

7 CONTINUE

C ПЕРЕСЫЛКА ЗНАЧЕНИЙ ИНТЕГРИРУЕМЫХ ПЕРЕМЕННЫХ В МАССИВ Y0(N)

DO 8 I=1,N

8 Y0(I)=Y(I)

11 RETURN

END

C ПОДПРОГРАММА ВЫЧИСЛЕНИЯ ПРАВЫХ ЧАСТЕЙ СИСТЕМЫ

С  3-ГО ПОРЯДКА:

C dy1/dt=f1(y1,y2,y3,x); где: f1(y1,y2,y3,x)=cos(x)-exp(-x)-y3

C dy2/dt=f2(y1,y2,y3,x);      f2(y1,y2,y3,x)=y3

C dy3/dt=f3(y1,y2,y3,x);      f3(y1,y2,y3,x)=-y2

SUBROUTINE PRAV(X,Y,F)

DIMENSION Y(*),F(*)

F(1)=COS(X)-EXP(-X)-Y(3)

F(2)=Y(3)

F(3)=-Y(2)

RETURN 

END