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

С       переменной x, которые пропускаются при выводе решения на экран

C      или печать (входная величина).

C      Например, при X0=0, XKAN=0.1, H=0.0001 решение

C      системы уравнений будут выссчитываться  в (XKAN-X0)/H=(0.1С      0.0)/0.0001=1000 узлах оси переменной x.

C      Но для анализа переходного процесса исследователь считает

C      достаточным сделать вывод решения системы  в 100

C      узлах. Тогда INTWR=1000/100=10 (будет  сделан вывод результатов

С      расчета в каждом десятом узле независимой переменной x);

C X - независимая переменная x (рабочая и выходная величина)

C X0 - начальное значение независимой переменной x (входная величина);

C XKAN – окончательное значение независимой переменной x (входная

С  величина);

C H - шаг изменения независимой переменной x (входная величина);

C Y(16) - массив текущих значений вектора Y(x) с компонентами  y(1),y(2),

C     ... y(j), ... ,y(n) (рабочие и выходные величины);

C Y0(16) - массив начальных значений вектора Y(x) (входные величины);

C Z1,Z2,Z3 - символьные переменные – имена файлов входных и выходных

C      данных

C NREZ – счетчик количества записей, которые выводятся в выходные   

С файлы

C PRAV - имя файла, в котором находится подпрограмма вычисления

С  правых частей

C        f(x,Y) системы дифференциальных уравнений. Это имя включено в

C        список фактических параметров оператора CALL RNGKT4(...,PRAV)       

С         этой программы. В списке формальных параметров 

C        оператора SUBROUTINE RNGKT(...,FNC) идентификатору  PRAV

С        соответствует второе имя FNC, то имя PRAV описано в операторе

C        EXTERNAL

DIMENSION Y(20),Y0(20)

EXTERNAL PRAV

CHARACTER Z1*8,Z2*8,Z3*8

C        Чтение имени файла входных данных с экрана дисплея

WRITE(5,*)'DIFEIL: ВВЕДИТЕ  ИМЯ ФАЙЛА ВХОДНЫХ ДАННЫХ',

*' В ВИДЕ  AAAA.DAT'

READ 1,Z1

1 FORMAT(A8)

WRITE(5,*)'DIFEIL: ИМЯ ФАЙЛА ВХОДНЫХ ДАННЫХ ВВЕДЕНО'

C        Открытие файла входных данных на диске и чтение входных

C        данных в оперативную память

OPEN(UNIT=1,FILE=Z1,STATUS='OLD')

READ(1,*)N,INTWR,X0,XKAN,H

READ(1,*)(Y0(I),I=1,N)

CLOSE(UNIT=1)

WRITE(5,*)'DIFEIL: ВХОДНЫЕ ДАННЫЕ ПРОЧИТАННЫЙ Из ДИСКА'

C        Образование имени двух файлов результатов в виде AAAA.REZ и

C        AAAA.GRF, где AAAA - имя файла входных данных.

С       В файл AAAA.REZ

C        программа DIFEIL засылает входные данные и таблицу решений

C        системы, которые обеспечены "шапками".

С        В файл AAAA.GRF засылается только таблица Y(N)

C        (выходная информация в таком виде нужна для программы

C        GRAF.EXE, которая преобразует  табличную форму

C        представления информации в форму осциллограмм на экране

С        дисплея)

Z2(1:4)=Z1(1:4)    

Z3(1:4)=Z1(1:4)

Z2(5:8)='.REZ'

Z3(5:8)='.GRF'

C        Открытие файлов результатов на диске

OPEN(UNIT=3,FILE=Z2,STATUS='NEW')

OPEN(UNIT=4,FILE=Z3,STATUS='NEW')

C        Запись входных данных в выходной файл AAAA.REZ

WRITE(3,*)'DIFEIL - РЕШЕНИЕ СИСТЕМЫ ДИФЕРЕНЦИАЛЬНЫХ',

*' УРАВНЕНИЙ МЕТОДОМ ЭЙЛЕРА'

WRITE(3,110)N,INTWR,X0,XKAN,H

110 FORMAT(1X,'ВХОДНЫЕ ДАННЫЕ'/1X,'         N     INTWR        X0

*      XK         H'/1X,2I10,3F10.5)

WRITE(3,111)(Y0(I),I=1,N)

111 FORMAT(1X,'Y0(N)'/(8G10.5))

C        Подготовка к реализации алгоритма методом Эйлера

X=X0                      

NREZ=0

C       Пересылка начальных значений Y0 интегрируемых переменных в

С       массив Y

DO 3 I=1,N

3 Y(I)=Y0(I)

C       Выдача заголовка таблицы результатов в выходной файл AAAA.REZ

WRITE(3,*)'           РЕЗУЛЬТАТЫ РАСЧЕТА'

WRITE(3,*)'      x    y(1)    y(2)    y(3)    y(4)    y(5)    y(6)

*    y(7)    y(8)'

IF(N.LE.8)GO TO 4

WRITE(3,*)'           y(9)   y(10)   y(11)   y(12)   y(13)   y(14)

*   y(15)   y(16)'

C Начало цыклической части решения системы уравнений

C        Надо ли записывать решения текущего шага в выходные файлы

C        AAAA.REZ i AAAA.GRF?

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

C        Увеличение счетчика записей на единицу

NREZ=NREZ+1