|
Главная -> Применение эвм 0 1 [2] 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 оператором DIMENSION А (4, 4), а фактический размер равен в лай ной задаче трем. Тогда девять элементов матрицы будут размещець в 16 зарезервированных под массив А последовательных ячейках памяти следующим образом; Ои «21 1-Й столбец 2-й Столбец 3-й столбец 4-й столбец Пример показывает, что если фактический размер матрицы ..VI меньше указанного макси.мального размера МО в операторе DIMENSION, то элементы в ячейках памяти располагаются не подряд. Порядковый номер ячейки К. соответствующей элементу А (I, .J), определяется формулой К = (-1 - I) * МО + I. где МО - максимально возможная длина столбца. Однако, как отмечалось выше, в вьоываемой стандартной подпрограмме матрица Л представлена в виде массива переменной длины, а элементы матрицы должны быть расположены в массиве А подряд без пропусков ячеек. Поэтому при обращении к стандартной подпрограмме, проводимом при М = 3, будут использованы числа, содержащиеся в первых девяти (М * М = 3 3 = 9) ячейках, зарезервированных под массив А в вызывающей программе. Очевидно, что эти числа не соответствуют коэффициентам построенной в вызывающей программе матрицы А (3, 3). Это происходит потому, что при описании матрицы в подпрограмме в виде массива А (Д\, М) переменной длины выбор номера К ячейки памяти, соответствующей элементу А (I, J), производится иа основе формулы К ~- (J - I) * М 4- I, где М - фактическая длина столбца, указанная при обращении к подпрограмме. Таким образом, при М ?ь д\0 матрица А будет передаваться в подпрограмму неправильно. Чтобы избежать рассмотренного несоответствия, следует перейти к векторной форме хранения матриц, при которой каждый столбец располагается в памяти непосредственно за предыдущим. В этом случае все М * М элементов расположатся в памяти подряд, а последующие (МО * МО - М * М) ячеек из общего числа зарезервированных останутся пустыми. Для нашего примера элементы матрицы А при М - 3 будут при использовании векторной формы хранения расположены так:
i-й столбец 2-й столбец 3-й столбец свободная область Переход к векторной форме хранения матриц, необходимый для работы со стандартными подпрограммами, можно осуществить одним из следующих трех способов. Можно сформировать матрицу в вызывающей программе в оной форме, а перед обращением к стандартной подпрограмме vniecTBHTb преобразование двумерной формы в векторную, т. е. I зовать «сжатие» области памяти, в которой находятся элементы Р Эта операция может быть выполнена с помощью специаль-ойстандартной подпрограммы ARRAY, которая осуществляет пре-обоазование одномерного массива матрицы в векторной форме в двумерный и наоборот. Обращение к этой подпрограмме имеет вид CALL ARRAY (MODE, N, М, NO, МО, В, А), где MODE - код преобразования (I - одномерный массив в двумерный, 2 - двумерный в одномерный); N. М - фактическое число строк и столбцов матрицы; N0, МО - число строк и столбцов, определенных для матрицы в операторе DIMENSION; В - одномерный массив длиной N * М (вектор); А - двумерный массив (матрица). Вектор В и матрица А могут быть в памяти совмещены. Например, рассмотренное выше преобразование матрицы А (4, 4) в векторную форму при М = 3 реализуется оператором CALL ARRAY (2, 3, 3, 4. 4, А, А). 2, Д10ЖИ0 формировать матрицу в подпрограмме, входными параметрами которой являются массив А, число строк N и число столбцов М. При этом в подпрограмме следует описать массив А как двумерный массив переменного размера с помощью оператора DIMENSION А (N, М). Тогда при формировании матрицы все ее элементы будут записаны подряд без «пробелов» в памяти, как это н требует векторная форма. В головной вызываюшей программе надо лишь определить максимально возможную длину массива А, причем его можно омнсать как одномерный, например DIMENSION А (1000). Описанный путь реализован в большинстве приводимых в дайной книге программ. 3. Можно сразу формировать матрицу в виде одномерного массива, в котором индекс К. соответствующий индексам I, J матрицы с N строками и М столбцами вычисляется по формуле К (J - I) * N + Если матрица имеет специальный вид. например является сим--четричной или ленточной, то для работы со стандартными подпрограм-ами коэффициенты матрицы всегда должны быть записаны подряд в одномерный массив в последовательности, зависящей от вида матри-ДЬ( и используемой подпрограммы. Примеры таких способов записи матриц будут рассмотрены ниже. Перейдем к краткому описанию ряда наиболее важных подпро-1"Рамм. Характеристика стандартных программ. Подпрограмма обращения пгрицы MINV реализует вычисление обратной матрицы А- ме- годом Гаусса-Жордана [2. Вызов этой подпрограммы производится так: CALL MINV(A, N, D, L, M), где A - входная матрица размером N X N, которая иа выходе заменяется обратной матрицей; D - определитель матрицы, вычисляемый в подпрограмме; L и М - целые рабочие массивы длиной N, которые следует описать в вызывающей программе. Подпрограмма решения системы линейных уравнений обиего вида GELG реализует решение методо.ч последовательного исключения Гаусса с выбором главного элемента. Эта и последующие подпрограммы предусматривают возможность решения N систем с одной и той же матрицей А, нос различными столбцами правых частей В. Для этого правые части задаются как матрица размером М / N. а N векторов решений также расположены в одном массиве последовательно по М элементов. Такая возможность реализована с целью экономии машинного времени, поскольку в случае N отдельных обращений к подпрограмме с разными правыми частями В над матрицей А будут (гроизводиться одни и те же операции исключения неизвестных. Обра-пеиие к подпрограмме имеет вид CALL GELG (В. А, М, N, EPS, lER), где В - матрица правых частей размером М х N, которая на выходе заменяется .матрицей решений; А - матрица коэффициентов системы уравнений размером М X М (в ходе вычислений эта матрица «портится»); М - число уравнений; N - число правых частей; EPS - входная константа, используемая при контроле потери точности; IER-выходной параметр ошибки. В программе в процессе реализации процедуры исключения предусмотрен контроль потери точности. Для этого на каждом шаге k исключения главный элемент аи, на который производится деление в формулах (1,11), (L14), сравнивается с величиной EPS * аях. тде "max - макси.чальный но модулю элемент матрицы. Если абсолютная величина главного элемента оказывается меньше, чем это произведение, то параметр ошибки принимает значение {к + !) Обычно (адают EPS = (10~ -f- Ю-**). Выходной параметр ошибки IER = " i в случае, когда на каком-то шаге исключения главный эле-MCfiT равен нулю и, следовательно, решение получено быть не может. Значение IER = О свидетельствует об отсутствии ошибок. Параметр IER следует анализировать после обращения к подпрограмме GELG, и в случае IER О предусматривать соответствующее сообщение. Для систем линейных уравнений с матрицами специального вида разработаны программы, реализующие модификации метода исключения с учетом структуры матрицы. При этом экономится память при хранении матриц, а в алгоритме исключения не проводятся опе- рации над нулевыми коэффициентами, что сокращает время решения системы. Программа GELS реализует решение системы линейных урав-пений с симметричной матрицей коэффициентов. Обращение к ней имеет вид CALL GELS (В, А, М, N, EPS, 1ER, AUX). Здесь А - массив, в который последовательно по столонам загш-саиа верхняя треугольная часть симметричной матрицы, л. е. этот массив содержит М * (М 1)/2 элементов в такой последователь- ности: Йп. 121 22, fljAf, .... ймм\ AUX - рабочий О о Ог О Рис !.3 массив длиной (М - 1), который должен быть описан в вызывающей программе; остальные параметры по смыслу полностью совпадают с одноименными параметрами подпрограммы GELG. При формировании массива А коэффициент матрицы (при i < /) записывается в элемент массива А (К), но.чер которого определяется по формуле К = J * - - 1)/2+ L Решение методом Гаусса систем линейных уравнений с ленточной матрицей реализует [юдпрограм.ча GELB. Существенной особенностью этой подпрограммы является способ записи коэффициентов матрицы. Для его пояснения рассмотрим матрицу с числом верхних диагоналей 2 и числом нижних диагоналей 1, показанную на рис. [.3. В массив А следует записать/70 ст/?ол:ал( все коэффициенты, лежащие в пределах ленты матрицы, т. е. для данного при.мера массив А должен содержать такую последовательность коэффициентов: йц, ai, й,з, Й22. 22., 0,3, Й24, Й54, Й55. Подпрограмма вызывается с помощью следующего обращения: CALL GELB (В, А, М, N, MUD, MLD, EPS, IER), гдеМ ~~ число уравнений; N - число правых частей; MUD верхних диагоналей, т. е. лежащих выше главной; MLD нижних диагоналей, т. е. лежащих ниже главной; А - матрица, записанная по строкам в одномерный массив, под который в вызывающей Программе должно быть отведено (М * МС ~ ML * (ML + 1)/2) 1; здесь ML = (МС - MLD - 1), МС = min (М, 1 + MUD J" MLD) - ширина ленты, остальные параметры совпадают с одно- енными параметрами подпрограммы GELG. gycVn 3 н 4 будут рассмотрены еще две подпрограммы: МГНр ~ решения систем с трехдиагональиой матрицей и ° ~-для решения систем уравнений с ленточной симметричной Р"ией методом квадратного корня. число число § 1.4. ПРОГРАММНАЯ РЕАПИЗАЦИЯ РАСЧЕТА СТАЦИОНАРНЫХ СРЕДНИХ ТЕМПЕРАТУР При составлении программы для расчета стационарных средних температур в системе тел н потоков теплоносителей используем систему уравнений (1.2)-(1.5). Приравнивая производные по времени нулю и подстаатяя выражения (1.5) для средней температуры теплоносителя Vi в уравнения баланса (1-2), (КЗ), получим систему алгебраических уравнений относительно температур Ti (i 1..... N.,) UT- Т" {I = 1, л/ж)- ДЛя составления программы целесообразно сгруппировать коэффициенты при неизвестных и записать эту систему в следующей форме: уравнения для твердых тел ii ~ 1, Nj) ,, /= 1 (=1 j /= 1 i-1 уравнения для объемов с теплоносителем (/ = 1, .... N 1=1 \ - I ЦВЫК (Ь27) уравнения для определения температур теплоносителя и а вхо-де Ul с \ Ж 2 c,,G„„+2c,G,, 2 c.n0iVT-- 2 0,G,,e,;(].28) здесь в уравиеиии (1.27) С, 2 G/ + 2G,t,. Представим систему (1.26)-{1.28) в матричной форме где А - матрица размером N х N, yV = yV -- 2 * Л,; X - столбец неизвестных температур, содержащий Ы членов; В - столбец из jV свободных членов. Примем следующий порядок нумерации jV неизвестных температур: первые членов столбца X соответствуют температурам тел Л.....Тл.. далее парами следуют 2 * N. температур теплоносителей на выходе и на входе Of, Будем использо- системы: первые yV, далее следуют парами также следующую нумерацию уравиеиии tuPHHfi совпадают с уравнениями (1.26). да RTypaBuenm 0.27), 0-28). Приводимая инже программа расчета средних температур , 1 4) оформлена в виде подпрограммы, параметрами которой вляются все необходимые исходные данные н рабочие массивы, писанные в комментариях к тексту. На основе принятой иу.мераи.ии лясистемы(1.26)-(1.28) в подпрограмме SYSTT произюдится формирование массивов А и Т, соответствующих матрице и столбцу свободных членов, а затем с помощью стаидартион подпрограммы GELG {см. § 1-3) решается система уравнений и находится столбец неизвестных. Отметим, что после решения системы найденные температуры размещаются в том же массиве Т, в котором находились свободные члены. В уравнениях баланса в общем случае могут быть связаны температуры Г,-, Tj, Ui с любыми ио.мерами /, /, /, поэтому матрица А является матрицей общего вида и для решения системы используется подпрограмма GELG. Рассмотрим процедуру формирования матрицы А н столбца Т. Сначала двумерный массив А и одномерный массив Т обнуляются, а затем производится расчет их ненулевых элементов путем последовательного суммирования отдельных членов, входящих в формулы (1.26)-(1.28). Организация этой процедуры суммирования зависит от используемого способа описания теплового взаимодействия между элементами системы. Для примера поясним способ описания теплового взаимодействия между телами. В случае, когда каждое тело взаимодействует со всеми остальными, все а]} отличны от нуля, и для описания тепловых связей можно использовать двумерный массив, элементами которого являются о]}. Однако обычно на практике число взаимодействующих пар тел значительно меньше максимально возможного значения, и поэтому большинство тепловых проводимостей равно нулю. С целью экономии машинной памяти и сокращения объема исходных данных целесообразно приводить информацию только об отличных от нуля Oi/- В приведенной программе эта информация задается с помощью двух одномерных массивов IJ и S1J. В массиве U длиной 2 * N1J, где NIJ - число связей между телами, парами записаны номера взаимодействующих тел: элементы [J {2 * М - I), [J {2 * М) указывают номера тел, участвующих в М-м взаи.модействни. Массив SIJ длиной NIJ содержит значения тепловых проводимостей а", соответствующих этим тепловым связям. Совершенно аналогичным образом представлена информация о тепловых связях ajf между телами и теплоносителями с помощью массивов IL и SIL, а также о связях gJI между телами и средами с ло.чощью массивов (К и SIK- Движение потоков теплоносителей описывается с помощью массива ML длиной 2 * NG и массива GC длиной NG, где NG - число по- 0 1 [2] 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 0.0135 |
|