ФЕДЕРАЛЬНОЕ АГЕНТСТВО ПО ОБРАЗОВАНИЮ

ГОСУДАРСТВЕННОЕ ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ

ВЫСШЕГО ПРОФЕССИОНАЛЬНОГО ОБРАЗОВАНИЯ

КЕМЕРОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ

Математический факультет

Кафедра Вычислительной математики

ДИПЛОМНАЯ РАБОТА

Разработка программного модуля для нахождения оптимальных предельно-допустимых выбросов в атмосферу от группы источников

студента 5 курса

Гурова Павла Владимировича

Кемерово 2009

ОГЛАВЛЕНИЕ

ОГЛАВЛЕНИЕ

ВВЕДЕНИЕ

1. Существующая система установления ПДВ для промышленных источников

1.1 Нормативы и показатели загрязнения атмосферы.

1.2 Данные об источниках загрязнения атмосферы.

1.3 Метеопараметры

1.4 Данные наблюдений за загрязнением атмосферы

1.5 Модели расчета загрязнения атмосферы

1.6 Унифицированные программы расчета загрязнения атмосферы

1.7 МетодЫ равного квотирования и МРН-87

2. Симплекс-метод

2.1 Общая характеристика симплекс метода

2.2 Алгоритм симплекс метода (первая симплекс таблица)

3. Формализация поставленной задачи

4. Программная реализация и пример практического применения

4.1 Выбор загрязняющих веществ

4.2 Обработка точек с повышенным загрязнением

4.3 Обработка источников

4.4 Обработка таблиц влияния источников на точки

4.5 Применении симплекс-метода

4.6 Вывод полученных результатов

4.7 Сравнение различных методов расчета ПДВ для реального предприятия

Заключение

Список литературы

Вспомогательные указатели

Перечень сокращений

Приложения

Unit1.pas

Simplex.pas

ВВЕДЕНИЕ

В нашей стране существует государственная система управления выбросами промышленных источников загрязнения атмосферы (ИЗА) [1-4]. Целью системы является достижение и поддержание на нормативном уровне N загрязнения воздушного бассейна населенных мест, которое характеризуется приземной разовой (осредненной за 20 минут) концентрацией С.

Основными разделами данной системы является проведение воздухоохранной экспертизы вновь строящихся объектов и установление предельно допустимых норм (ПДВ) для существующих объектов.

В этих задачах расчетным путем устанавливается связь между выбросами Q предприятий и приземными концентрациями C загрязняющих веществ (ЗВ), выбрасываемыми промышленными источниками.

Для расчета концентраций С по известным значениям Q и характерным для данной территории метеопараметрам используется единая, утвержденная на государственном уровне методика ОНД-86 [2].

В настоящее время проектные расчеты осуществляются с использованием программных комплексов (ПК) сертифицированных главной геофизической обсерваторией им. А.И. Войекова (г. Санкт-Петербург).

При проведении атмосферноохранной экспертизы проекта нового предприятия, решается прямая задача моделирования загрязнения атмосферы.

При установлении ПДВ для существующего предприятия, решается обратная задача – подбираются такие значения выбросов Х, при которых в прилегающих к предприятию жилых районах приземные концентрации С не превышают санитарную норму N.

В соответствии с нормативными документами по нахождения ПДВ [1,4] в этом случае обратная задача решается путем подбора реальных (технологически применимых) для данного производства атмосфероохранных мероприятий (АОМ), таких как замена оборудования, установка газоочистки и т.д. Фактически, обратная задача решается при этом методом подбора. Однако при большом числе источников большую помощь в выборе первоочередных мероприятий дает решение задачи форматизированными методами. Пользуясь линейностью модели ОНД-86 по выбросам источников можно представить загрязнение атмосферы в контрольных точках жилой зоны в виде линейной формы:

,

где - концентрация в i-ой точке, - выброс j-го источника, - вклад j-го источника в i-й точке, который в дальнейшем будем называть коэффициентом влияния.

Санитарные требования приводят к системе линейных неравенств:

,

решение которой ищется на интервале

,

где - технологически обоснованный минимум выброса j-го источника.

На сегодняшний день в методической литературе описаны два метода нахождения решения поставленной задачи: МРН-87 [24] и метод равного квотирования [25]. Оба метода дают частное решение поставленной системы неравенств из соображений удобства нахождения единственного решения. Однако любое предприятие заинтересовано в минимальных затратах, необходимых для установления нормативных выбросов . Для этой цели, к поставленной системе неравенств добавляем целевую функцию:

,

где в общем случае - стоимость снижения на единицу выброса для j-го источника. В данном виде решение X_j дает минимум затрат на достижение нормативного загрязнения атмосферы.

В предположении получаем

,

что эквивалентно .

Следовательно, задача сводится к поиску максимально допустимого по сумме сочетания выбросов X_j данного предприятия, позволяющего достичь нормативного загрязнения атмосферы. В связи с линейностью модели ОНД-86 по отношению к выбросам для поиска решения поставленной задачи может быть применён симплекс-метод.

Целью настоящей работы является разработка программного модуля для расчета оптимальных значений ПДВ с использованием симплекс-метода. Данный модуль должен не просто реализовать решение задачи линейного программирования, но и быть интегрирован в среду программного комплекса ЭРА-ВОЗДУХ, который предназначен для решения проектных задач, использует стандартизированные данные об источниках выбросов, метеорологических параметрах распространения примесей и выдает выходную информацию в установленной нормативными документами форме (подробности на сайте www.logos-plus.ru). После разработки и отладки модуля в работе предполагается решить задачу на основе реальных данных и сравнить оптимальное решение с методами, описанными в методической литературе.

1. Существующая система установления ПДВ для промышленных источников

1.1 Нормативы и показатели загрязнения атмосферы

Основными нормативами чистоты атмосферы, лимитирующими загрязнение (следовательно и выбросы ИЗА), являются предельно допустимые концентрации (ПДК) для населенных мест, которые подразделяются на среднесуточные (ПДКс - длительных периодов воздействия) и максимально-разовые (ПДКр - 20-ти минутного осреднения) [6,7,8]. Такие ПДК отражают критерии качества воздуха для человека (гигиенические). В документах [1,5] при определении допустимых выбросов в атмосферу предусмотрен учет ПДК для растительности, животного мира, сельхозугодий и т.д. (экологические) с выбором в качестве определяющей наиболее жесткой величины. Однако пока утвержденные какими-либо нормативными документами экологические ПДК отсутствуют [1], а управление выбросами ведется с ориентацией только на гигиенические ПДК.

ПДКс имеет цель предупреждение хронического воздействия ЗВ при длительном вдыхании. Концентрация ЗВ в пределах ПДКс не должна оказывать на человека прямого или косвенного воздействия при неограниченно продолжительном вдыхании. Поэтому и гигиенистами [9] и специалистами по атмосферной диффузии [10] отмечается, что ПДКс, несмотря на название, можно использовать для таких времен осреднения, как год или сезон. Список ПДКс дополняют максимально разовые предельно допустимые концентрации (ПДКр), при установлении которых учитывается рефлекторное воздействие. В случае отсутствия какого-либо ПДК используется уравнение ПДКр=10ПДКс. Для ряда ЗВ, по которым пока не определены нормативы ПДК, расчетным путем установлены ориентировочно безопасные уровни воздействия (ОБУВ), используемые в качестве ПДКр. На сегодняшний день в перечне [7] около 3000 веществ имеющих какой-либо норматив ПДК. В нашей стране в качестве основной размерности ПДК принят мг/м³ .

Рассмотрение списка ПДК показывает, что отношение ПДКр/ПДКс лежит в интервале от 1 до 20.

Наличие ПДК, регламентирующих как длительное, так и кратковременное воздействие ЗВ на здоровье человека, свидетельствует, что показателями загрязнения атмосферы населенных мест должны являться концентрации, согласованные с ПДКс и ПДКр по времени осреднения. Таким образом, названными показателями будут среднегодовая Сг и разовая С концентрации, периодами осреднения которых являются год и 20 минут соответственно.

В настоящее время, согласно [1], считается, что чистота ВБ населенных мест для некоторого ЗВ соответствует нормативной, если на границах санитарно-защитных зон предприятий и в пределах жилой зоны выполняется неравенство

См < mПДКр (1.1)

где См – максимальная суммарная концентрация данного ЗВ; m = 0.8 для особо охраняемых территорий и m = 1.0 в остальных случаях. К особо охраняемым территориям относятся места массового отдыха населения (подробно см. [16]) и территории размещения лечебно-профилактических учреждений длительного пребывания больных и центров реабилитации. В качестве См берется вычисленная по [2] суммарная от всех источников концентрация, достигающаяся при неблагоприятном сочетании метеорологических условий (МУ) и параметрах выбросов ИЗА данного ЗВ. Однако, как отмечается в [10], “при планировании долгосрочных АОМ часто невозможно или экономически нецелесообразно предусматривать, чтобы абсолютно ни при каких МУ концентрации примесей не превышали ПДКр. Обычно оказывается достаточным учесть часто наблюдаемые неблагоприятные МУ, а на периоды особо опасных аномальных МУ вводить в действие оперативные АОМ”. Поэтому в неравенстве (1.1) под С_м подразумевается не абсолютный максимум, а разовая концентрация при так называемых “нормальных неблагоприятных метеоусловиях”. Частота появления МУ, способных вызвать большее загрязнение ВБ (называемых “аномальными, неблагоприятными”) при тех же уровнях выбросов зависит от метеорологического режима района расположения предприятий [10]

Отсутствие в нормативных указаниях [1] аналогичного (1.1) неравенства для годового осреднения

Сг < ПДКс (1.2)

объясняется тем, что если (1.1) выполняется, то (1.2) для того же ЗВ будет следовать автоматически. Действительно, отношение ПДКр/ПДКс за редким исключением не превосходит 10, в то врем как отношение См/Сг числом из интервала от 25 до 40 [10]. Однако этот интервал относится прежде всего к случаю изолированного или группы близко расположенных ИЗА. В крупном промышленном городе, при наличии большого числа распределенных по территории и различных по типу источников, нижняя граница указанного интервала существенно снижается. Оценки отношения См/Сг на основе статистического анализа результатов стационарных наблюдений в различных городах России рассматриваются в работах Э.Ю.Безуглой [11,12]. Фактические соотношения между См/Сг для серистого ангидрида, полученные по данным регулярных наблюдений на 863 стационарных постах, расположенных в городах России показывают, что нормирование выбросов по разовым концентрациям в 30% случаев не обеспечит соблюдение (1.2). Кроме того, следует подчеркнуть, что определение риска от загрязнения атмосферы для здоровья населения основано именно на долгопериодных концентрациях ЗВ. Отсюда следует, что в общем случае следует предусматривать проверку обеих условий (1.1) и (1.2).

С начала 90-х годов в состав нормативных данных по загрязняющим веществам добавились базовые удельные стоимостные показатели [19] по выбросам каждого ЗВ за превышение нормативов предельно допустимых выбросов (ПДВ). Под ПДВ для предприятия понимается такая совокупность выбросов его ИЗА [1,3,4,5], при которой (с учетом фонового загрязнения атмосферы от других предприятий) соблюдается условие (1.1). За выброс в пределах установленного проектом нормативов ПДВ лимита предприятие платит базовый норматив (руб/тонну в год), а за превышение – в пятикратном размере. Кроме того, начисляемая плата умножается на коэффициент экологической ситуации (территориальный показатель, от 1.0 до 1.44) и коэффициент индексации. Для примера можно сказать, что годовая плата бытовой котельной может составлять от нескольких тысяч до нескольких десятков тысяч рублей, а крупного химического или металлургического предприятия достигает нескольких миллионов.

Здесь хотелось бы еще раз подчеркнуть, что при разработке проекта норм ПДВ условие (1.1) проверяется расчетным путем на основе данных о выбросах ИЗА с использованием методики ОНД-86 [2]. Следовательно, и плата за выброс зависит от локальной модели расчета загрязнения атмосферы, которая нормативно утверждена на федеральном уровне. Как и любая (а тем более - инженерная) модель столь сложного явления как атмосферная диффузия, она не описывает всех природных закономерностей и может допускать значительные погрешности. Тем не менее данная проблема в нормативных работах не обсуждается, поскольку в разрешении финансовых вопросов обязательно требуется ссылка на законодательно утвержденный документ, которым и является ОНД-86.

1.2 Данные об источниках загрязнения атмосферы

Естественно, важнейшей информацией для управления выбросами промышленных предприятий являются параметры их ИЗА. Существует ориентированная на данные цели система сбора и хранения такой информации на государственном уровне. При постановке многих научных и экспертных задач, не входящих в перечень нормативных проектных проработок, исследователями явно или неявно предполагается, что названная система обеспечивает (или обеспечит в ближайшем будущем) представление информации об ИЗА в необходимом пространственном и временном разрешении и с требуемой для поставленной задачи точностью. В этой связи представляется интересным рассмотреть информацию об ИЗА и те ограничения которые должны учитываться не только при постановке, но и при решения исследовательских задач применительно к реальным объектам промышленных городов или регионов.

Источник загрязнения атмосферы (ИЗА) представляет собой хозяйственное сооружение или естественный объект, из которого в атмосферу поступают ЗВ. К источнику загрязнения, как правило, относится один или несколько источников выделения, которыми называются промышленные установки, генерирующие ЗВ в процессе работы (например, водогрейный котел есть источник выделения, а труба – источник загрязнения).

Нормативная информация о параметрах ИЗА представляется предприятиями в природоохранные органы на основе проведения инвентаризации. В Федеральный Закон "Об охране атмосферного воздуха" [5] впервые (в сравнении с ранее действовавшим Законом) введена статья, касающаяся инвентаризации выбросов ЗВ в атмосферный воздух. В соответствии со статьей 22 Закона “юридические лица, имеющие источники выбросов ЗВ в атмосферный воздух, проводят инвентаризацию выбросов вредных (загрязняющих) веществ в атмосферный воздух и их источников в порядке, определенном специально уполномоченным федеральным органом исполнительной власти в области охраны атмосферного воздуха» (в настоящее время - Ростехнадзор РФ). Инвентаризацию проводят все действующие предприятия, организации, учреждения независимо от их организационно-правовых форм и форм собственности, производственная деятельность которых связана с выбросом ЗВ в атмосферу. Ответственность за полноту и достоверность данных инвентаризации несет предприятие (в лице руководителя). Инвентаризация выбросов (т.е. представление информации с точностью до каждого ИЗА) проводится 1 раз в 5 лет [1, 15]. Кроме того, ежегодно предприятия отчитываются по форме государственного статистического наблюдения № 2-тп (воздух), которая содержит сведения только о годовых выбросах предприятия в сумме по всем источникам.

Определение параметров источников загрязнения атмосферы (ИЗА) для проекта нормативов ПДВ должно осуществляться при регламентных загрузке и условиях эксплуатации технологического и пыле-, газоочистного оборудования. В соответствие с новыми требованиями [1], параметры ИЗА следует фиксировать и на основных режимах работы технологического оборудования (установки) и различных стадиях технологических процессов. Однако в настоящее время для подавляющего большинства предприятий такая информация отсутствует, и данные о выбросах представлены только для максимальной регламентной нагрузки установок (максимальные разовые выбросы, г/сек) и в целом за год (годовые или валовые выбросы, т/год).

Для определения количественных и качественных характеристик выделений и выбросов загрязняющих веществ (3В) в атмосферу используются инструментальные и расчетные (балансовые, а также основанные на удельных технологических нормативах или закономерностях протекания физико-химических процессов) методы [1]. К расчетным методам, как правило, относятся расчетно-аналитические методы, в которых в качестве параметров расчетных формул для определения величин выброса (г/с) используются значения многократно измеренных концентраций вредных веществ (мг/м³ ) в атмосферном воздухе для типовых источников выделения. Выбор методов определения количественных и качественных характеристик выделений и выбросов 3В в атмосферу зависит, в первую очередь, oт характера производства и типа источника.

Инструментальные методы являются превалирующими для ИЗА с организованным выбросов ЗВ в атмосферу (организованные источники), к которым, в основном, относятся:

· дымовые и вентиляционные трубы, вентиляционные шахты;

· аэрационные фонари;

· дефлекторы.

Расчетные методы применяются, в основном, для определения характеристик источников с неорганизованными выделениями (выбросами). К ним относятся:

· неплотности оборудования (в т.ч. работающего при избыточном давлении);

· погрузочно-разгрузочные работы;

· открытое хранение сырья, материалов и отходов;

· оборудование и технологические процессы как в производственных помещениях, не оснащенных вентиляционными установками, так и расположенные на открытом воздухе (например, передвижные сварочные посты, резервуары хранения нефтепродуктов и т.д.);

· пруды-отстойники и накопители, нефтеловушки, шлако- и хвостохранилища, открытые поверхности испарения и т.п.; взрывные работы; открытые стоянки автотранспорта;

· передвижные источники, эксплуатируемые на производственной территории (автотранспорт, тепловозы, дорожная и строительная техника, речные суда и т.п.).

При этом, как подчеркивается в [1], могут использоваться только методики, рекомендованные к применению в установленном порядке перечнем [13] или документами, дополняющими и корректирующими этот перечень.

В настоящее время данные инвентаризации содержат сведения о типе источника, координатах (x1,y1), (x2,y2) [м] его расположения на плоской карте (схеме) местности, высоте выброса H [м], геометрии выходного канала D [м]. Кроме того, указывается выбрасываемый в единицу времени объем Vg [м³ /сек] (или вертикальная скорость Wo[м/сек]) газовоздушной смеси, ее температура Tg [^о С], показатели разового выброса Qi [г/сек] и суммарного за год (валового) выброса Qгi [т/год], i=1,…,I, где I - число выбрасываемых в атмосферу ЗВ.

Источники по типу подразделяются на точечные, линейные и площадные. Для точечного источника (x2,y2)=0. Линейные представляются прямолинейными отрезками конечной длины, при этом (x1,y1), (x2,y2) есть координаты концов отрезка. Площадные описываются прямоугольником с центром в точке (x1,y1) и длиной сторон (x2,y2). Для площадного ИЗА указывается острый угол, отсчитываемый против часовой стрелки от оси x до x2. В качестве показателей выбросов Qi и Qгi для линейных и площадных источников указываются суммарные значения, что предполагает равномерное выделение с единицы длины или площади. Источники произвольной формы приближаются набором отрезков или прямоугольников. Аналогично приближаются источники с неравномерным выбросом по протяженности или площади.

По величине погрешности данные об источниках можно разделить на две группы. Первая - характеристики источника как сооружения (высота H, диаметр D) и координаты его расположения на территории. Измерение таких данных осуществляется достаточно просто, в зависимости от потребностей задачи они могут быть легко уточнены и изменяются только при коренной перестройки промышленного предприятия. Поэтому при постановке задач можно предполагать, что перечисленные параметры ИЗА, взятые по результатам инвентаризации, заданы с малой погрешностью даже для локальных задач переноса аэрозолей. Вторая группа – характеристики условий и количества выхода ЗВ из источника: Tg, Vg (Wo), Q. Температура выброса является весьма устойчивой величиной для однотипных промышленных установок и погрешность ее задания на превосходит 2-3% для энергетических установок и 10-15% для остальных промышленных ИЗА. Большей изменчивостью обладает объем Vg (или скорость Wo газо-воздушной смеси), которые по имеющимся оценкам может изменяться для некоторых ИЗА в 2 раза. Такой важный параметр источника как выбрасываемое количество ЗВ в единицу времени Q также обладает значительной изменчивостью, и точность его получения весьма невысока. Реальные данные по конкретным предприятиям в большинстве случаев получают на основе расчетных методов, поэтому теоретически их погрешность неизвестна. На практике же специалисты оценивают ее величиной 50-100% для конкретных источников и 20% для предприятий в целом [20].

Из приведенного краткого рассмотрения следует, что нормативная информация о выбросах промышленных предприятий, представленная в инвентаризации, пригодна для моделирования максимального (причем без указания его наступления во времени) или среднего за достаточно длительный (сезон, год) период загрязнения атмосферы промышленными аэрозолями. Поскольку зависимость мощности разового выброса Q от времени не может быть в масштабе города (региона) напрямую получена по данным инвентаризации, моделирование временного хода загрязнения с учетом зависимости выбросов от времени может быть осуществлено лишь весьма приближенно, на основе косвенных данных по суточной или сезонной производительности промышленных объектов. Исключения составляют случаи специальных экспериментов, во время которых ИЗА оборудуются измерительной аппаратурой для контроля выбросов. А поскольку количество источников по многим ЗВ в промышленном городе измеряется сотнями, то очевидно, что такие эксперименты реальны только для уникальных ЗВ, допускающих надежные методы оперативного инструментального контроля и выбрасывающихся в атмосферу города незначительным числом ИЗА.

Можно также обоснованно полагать, что значения выбросов, полученные по данным инвентаризации в масштабе города (региона), если и содержат ошибку, то систематическую. Это связано с тем, что все выбросы получены на основе единообразных методов расчета, а случайные ошибки с определенной надежностью фильтруются контролирующими органами на этапе приемки и согласования инвентаризации. Тем не менее, даже при наличии систематической ошибки данные инвентаризации остаются пригодными для сравнительного анализа воздействия различных объектов на загрязнение атмосферы.

Что касается постановки и решения обратных задач по управлению выбросами в атмосферу, то тут следует учитывать как высокую чувствительность решения таких задач к погрешностям исходных данных, так и ограниченные возможности ИЗА реальных предприятий по оперативному изменению выбросов. Большинство промышленных установок не допускают отклонений от регламентного режима. Остановка многих из них либо вообще невозможна, либо сопровождается увеличением выбросов в атмосферу. Поэтому практическая ценность математических моделей (например - оптимизационных) оперативного управления разовыми выбросами весьма ограничена. А реальную пользу по планированию первоочередных АОМ могут принести модели оптимального планирования долгосрочных АОМ, где в качестве управляемых параметров являются суммарные годовые и максимальные за год разовые выбросы ИЗА.

1.3 Метеопараметры

Состояние локального загрязнения приземного слоя воздуха существенно зависит от метеорологических условий. Хорошо известно, что при одних и тех же параметрах выбросов ИЗА, в зависимости от метеоусловий, концентрация у земли может меняться на порядок и более.

С точки зрения распространения ЗВ в атмосфере метеоусловия подразделяются на нормальные и аномальные [10]. Нормальные характеризуются, прежде всего, наличием ярко выраженного среднего направления ветра. Таковыми в крупных городах являются условия со скоростью ветра более 1-2 м/сек. При меньших скоростях (штиль или близкое к штилевому состояние) в результате рельефных особенностей и температурной неоднородности подстилающей поверхности могут образовываться локальные циркуляционные зоны, приводящие к накоплению ЗВ в слое дыхания [17,18]. Ситуация становится особенно опасной при наличии вертикальной температурной инверсии, препятствующей уносу примеси в верхние слои атмосферы. Именно при таких метеоусловиях фиксируются максимальные уровни загрязнения при инструментальных наблюдениях.

Подразделение метеоусловий на нормальные и аномальные играет важную роль для осознания результатов инженерных расчетов загрязнения атмосферы. Дело в том, что все инженерные модели применимы только при нормальных метеоусловиях, поскольку единое направление ветра и его стационарность являются их непременным условием. Поэтому расчетная "максимальная" концентрация является не абсолютным максимумом загрязнения, а наибольшей из концентраций для нормальных условий. Даже если предположить, что методика расчета и параметры выбросов полностью соответствуют происходящему в природе, то превышение расчетного максимума все является равновозможным и зависит от частоты появления аномальных неблагоприятных метеоусловий.

В рамках любой локальной стационарной модели наиболее важными с точки зрения рассеяния примесей метеорологическими параметрами являются скорость и направление ветра, а также показатели диффузионной активности (устойчивости) атмосферы. Скорость и направление ветра измеряются непосредственно. Обзор методов измерений скорости ветра показывает, что относительная погрешность составляет от 15% (при скоростях порядка 5 м/с) до 55% (при скорости порядка 1м/с) [26]. С точки зрения решения задач переноса аэрозолей локального масштаба представляет интерес, что погрешности при определении направления ветра могут привести к ошибке положения оси дымового факела на карте территории порядка 10-15%, в результате чего при достаточно устойчивой стратификации атмосферы факел просто не накроет расчетную точку и приведет к большой ошибке моделирования. Это следует учитывать при интерпретации понятия "опасное направление ветра" и определении основных виновников загрязнения заданной точки города. Сказанное еще раз подчеркивает, что на практике при решении краткосрочных задач нормирования выбросов в смысле неравенства (1) представляет интерес предсказание с разумной точностью максимума разовой концентрации даже без указания момента времени, когда это произойдет.

1.4 Данные наблюдений за загрязнением атмосферы

При обосновании системы мониторинга на территории России академик Ю.А. Израэль подчеркивает, что "только регулярные наблюдения в строго определенных местах и в строго установленные сроки являются источниками прямой и статистически обеспеченной информации о загрязнении окружающей среды"[23]. Такого рода наблюдения применительно к загрязнению атмосферы осуществляются на сети стационарных постов Росгидромета в городах . Посты расположены только в крупных городах. Например, в Кемеровской области стационарные посты оборудованы в Новокузнецке (10), Кемерово (9) и Белово (1). Причем количество постов сокращается (в 90-х годах в Кемерово было 12 постов). На постах ежедневно (в 7, 13 и 19 часов местного времени) осуществляются отборы проб воздуха, которые доставляются в лабораторию местного подразделения Росгидромета, где анализируются стандартизированными методами. Контролируются только незначительное число наиболее распространенных ЗВ (10-20 примесей), в то время как в данных инвентаризации совокупности предприятий крупного промышленного города встречается на практике 100-200 веществ. Таким образом, большинство ЗВ на сети не контролируется ничем, кроме интуиции разработчиков проектов и согласующих эти проекты экспертов.

Обзор методов инструментального анализа воздуха [21,22] показывает, что количественные оценки погрешностей различных этапов лабораторных методов анализа полученных проб составляют 6-25% (с доверительной вероятностью 95%).

1.5 Модели расчета загрязнения атмосферы

Даже при наличие данных наблюдений за загрязнением воздуха на сети стационарных постов, одной из важнейших наукоемких задач охраны атмосферы является расчет загрязнения заданной территории по имеющимся данным о выбросах и условия распространения примесей. Действительно, оценить перспективный уровень загрязнения в зависимости от варианта промышленного развития можно только расчетными методами. Кроме того, методы инструментальных наблюдений в общем случае не могут указать вклад отдельного источника (предприятия) в измеряемую суммарную величину, что необходимо при определении основных виновников загрязнения, установлении ПДВ и начислении платы за выброс. Поэтому моделирование загрязнения атмосферы необходимо как для анализа, так и для прогноза.

В предисловии к [17] А.М. Яглом подчеркивает, что расчет диффузии примеси в атмосфере "не может быть сведен к какой-то задаче математической физики, а обязательно требует тех или иных нестрогих гипотез и приближенных допущений. По этой причине задача о распространении примеси в атмосфере не имеет одного общепринятого ‘правильного решения’, а характеризуется наличием целого ряда различных подходов к требуемому расчету, ни один из которых не может претендовать на полную строгость и точность".

Модели, естественно, можно классифицировать с различных точек зрения, и этому посвящена весьма обширная литература [10,17,18,20]. С позиций, используемых для построения научных теорий, модели подразделяются на статистические и полуэмпирические. Статистические основаны на том предположении, что поступающее из источника ЗВ переносится вместе со средним потоком, а его распространение в поперечном потоку направлении происходит под воздействием вихрей, движение которых подчиняется определенным статистическим закономерностям. Полуэмпирические используют для получения результата те или иные решения уравнения турбулентного переноса в предположение об аналогии между турбулентной и молекулярной диффузией. Слово "полуэмпирические" подчеркивает, что для задания коэффициентов турбулентной диффузии необходимы эмпирические предположения, и только после этого можно начинать поиск точного или численного решения уравнения переноса.

В зависимости от времени осреднения рассчитываемой концентрации локальные модели можно подразделить на краткосрочные (разовые) и долговременные. Как статистические, так и полуэмпирические краткосрочные модели рассчитывают концентрацию атмосферной примеси, осредненную за 20-30 минут. Долгосрочные предназначены для оценки загрязнения, осредненного за большие промежутки времени (сезон, год) и основаны на осреднении разовых расчетов с использованием повторяемости характерных для данной территории условий распространения ЗВ [17,24,27,31].

С точки зрения простоты использования модели также можно разбить на два класса: научно-исследовательские и инженерные. Первые, будучи способными описывать достаточно тонкие особенности распространения ЗВ в атмосфере (сложный рельеф, штиль, особенности турбулентного режима, локальные циркуляции и т.д.) [17,18], являются безусловно более сложными, требуют высокой квалификации персонала как для применения моделей, так и для получения специфических и дорогостоящих исходных данных. Вторые, предназначенные для проектных расчетов, доведены до однозначно трактуемых числовых зависимостей, табличных и графических аппроксимаций и обеспечены системой сбора (расчета) исходной информации для возможности реализации необходимых количественных оценок в процессе выполнения проектных работ [1,2]. Можно сказать, что именно наличие реальной на сегодняшний день системы обеспечения модели исходными данными является признаком того, что сама модель представляет практический интерес для целей управления (нормирования) промышленных выбросов.

Наиболее распространенной для инженерных приложений и принятия административных решений за рубежом является локальная модель Гауссовского факела [17,20], в основе которой лежит статистическая теория с рядом упрощающих предположений и большим количеством эмпирических таблиц для задания дисперсионных коэффициентов [20], определяющих процесс расширения факела при удалении от источника. Краткосрочная модель позволяет рассчитать в заданной точке территории разовую (среднюю за 20-30 минут) концентрацию в зависимости от трех параметров – скорости ветра, направления ветра и класса устойчивости (температурной стратификации) атмосферы. Гауссовская модель долговременного осреднения позволяет оценить среднюю за заданный период концентрацию по совместному распределению (за этот период) указанных трех параметров. Для практического использования модели достаточно иметь данные стандартных метеорологических наблюдений и параметров источников в объеме, рассмотренном в пункте 2.

Принятая на федеральном уровне Минприродой России [2] для нормативных расчетов модель атмосферного переноса ЗВ предназначена для расчета максимальных разовых концентраций См атмосферных примесей, фигурирующих в условии (1.1). Модель была создана в конце 60-х годов и в дальнейшем доведена до инженерных формул коллективом разработчиков под руководством проф. М.Е. Берлянда в Главной геофизической обсерватории им. А.И. Воейкова. В настоящее время она оформлена в виде руководящего документа ОНД-86 [2], вышедшего взамен ранее действующего СН-369-74 (с рядом уточнений, учетом коэффициента рельефа и застройки). Модель основана на численном решении стационарного полуэмпирического уравнения турбулентной диффузии для различных типов источников с последующей аналитической аппроксимацией. Численные расчеты проводились при атмосферных стратификациях, которым соответствуют наибольшие значения максимума приземной концентрации, достигаемые при “опасной” скорости ветра. Отсюда следует, что модель [2] изначально нацелена на расчет верхней оценки приземных концентраций при нормальных (соответствующих условиям стационарности) метеоусловиях. Поэтому аналитические аппроксимации, в отличие от Гауссовских моделей, позволяют оценивать концентрации ЗВ только в зависимости от двух параметров – скорости и направления ветра. А поскольку невозможно получить априорного аналитического решения для нахождения максимума (по скорости и направлению) суммарной концентрации для произвольного множества произвольно расположенных ИЗА, то указанная задача решается численно для каждого конкретного случая либо методом простого перебора, либо на основе анализа характерных особенностей заданного множества ИЗА. Для применения модели нужно знать, кроме параметров ИЗА, коэффициент рельефа территории, коэффициент осаждения исследуемой примеси, параметр характерной “опасной” стратификации атмосферы (все три этих параметра специально ориентированны на модель [2]) и климатические показатели, определяющие для заданной территории режим температуры окружающего воздуха и скорости ветра. Соответствующие количественные характеристики, вместе с фоновыми концентрациями, выдаются для проектных расчетов территориальными подразделениями Росгидромета.

1.6 Унифицированные программы расчета загрязнения атмосферы

С середины 70-х годов начали разрабатываться программы расчета загрязнения атмосферы (ПРЗА) для широко распространенных в то время ЭВМ типа БЭСМ и М-20. Наличие достаточно сложной вычислительной основы при расчете поля приземных концентраций и многовариантности методики [2] приводили к естественным ошибкам в алгоритмах и программных реализациях и, как следствие, к несопоставимости результатов, полученных на основе различных программ. Поэтому сложилась практика тестирования разработанных в различных организациях программ расчета экспертами ГГО им. А.И. Воейкова. Такая практика сохранилась и до настоящего времени. Успешно прошедшие тестирование программы получают статус "согласованных" и относятся к классу унифицированных ПРЗА (или УПРЗА). Список таких программ ежегодно рассылается Минприродой России своим территориальным подразделения и всем заинтересованным организациям. В настоящее время наибольший интерес представляют УПРЗА для персональных компьютеров, имеющие дружественный интерфейс и широкие графические возможности. Одна из УПРЗА входит в состав программного комплекса ЭРА-ВОЗДУХ (подробнее см. WWW.logos-plus.ru), разработчиком алгоритма и программной реализации которой является А.А. Быков.

Использование УПРЗА значительно облегчает расчет загрязнения приземного слоя атмосферы от существующего или проектируемого промышленного объекта. Для этого достаточно задать определяющие климатические параметры города или района расположения объекта, описать состояние источников загрязнения атмосферы (см. пункт 2.2), выбрать расчетный прямоугольник, задать параметры поиска максимума по скорости и направлению ветра. После этого надо запустить программу на расчет. В итоге пользователь получает в заданных расчетных точках таблицу максимальных концентраций для всех выбрасываемых объектом ЗВ и построенную на ее основе карту-схему предприятия с уровнями загрязнения атмосферы в виде изолиний. В таблице обязательно указываются опасные скорости и направления ветра, а также несколько, так называемых, "основных вкладчиков", на которые (с точки зрения модели) следует в первую очередь нацелить атмосфероохранные мероприятия. Кроме того, прилагаемые к УПРЗА сервисные блоки формируют всю необходимую проектную документацию, соответствующую по форме и содержанию единым требованиям.

Основным режимом любой УПРЗА, реализующих нормативную методику [2], является поиск в каждой заданной расчетной точке максимума приземной концентрации при всевозможных скоростях и направлениях ветра. Скорость варьируется в интервале от 0.5 до U*, где U* - скорость, вероятность превышения которой не более 5%. Поиск максимизирующего (или “опасного”) направления ветра осуществляется по всему кругу или в заданном пользователем секторе, если направление ветра, при котором предприятие влияет на жилые районы, достаточно очевидно. Естественно, можно рассчитать и поле приземных концентраций при заданных скорости и направлении ветра. Однако следует подчеркнуть, что использование такого расчета для целей оперативного мониторинга загрязнения не вполне корректно, поскольку стратификация в [2] априорно предполагается неблагоприятной, а не соответствующей данному моменту времени. Эту особенность следует учитывать и при сопоставлении более сложных исследовательских моделей с нормативной. Например, модели, основанные на достаточно реалистичных приближениях пограничного слоя атмосферы могут учитывать криволинейность воздушного потоков и продемонстрировать, что при ветре такого-то направления промышленный факел существенно смещается в сторону по сравнению с прямолинейным, полученным по расчету на УПРЗА при том же направлении. Безусловно, данный факт является основанием для серьезной критики модели, если ее целью является расчет “актуальных” концентраций для заданного момента времени. Но методика [2] и реализующие ее УПРЗА на ставят перед собой такой цели, а предназначены для получения верхней оценки разовой концентрации в каждой расчетной точке, полученной в процессе перебора возможный скоростей и направлений ветра. Поэтому основанием для критики [2] (цель не соответствует результату) со стороны разработчиков физически более совершенных моделей может служить расчет в той или иной области города более высоких концентраций при каких-то характерных условиях распространения ЗВ (естественно, при идентичных параметрах ИЗА).

1.7 МетодЫ равного квотирования и МРН-87

В соответствии с [43] определение допустимых выбросов X_j методом равного квотирования выполняются на основе анализа рассчитанных приземных концентраций C_i в контрольных точках, где расчетная модель ОНД-86 дает превышение ПДК.

Алгоритм нормирования

Рекомендуемое значение г/с выбросов рассчитываются исходя из принципа равного квотирования вкладов ИЗА, обеспечивающих концентрацию в точке в пределах ПДК. Квота вклада C_{н j} каждого j-го источника определяется для каждой точки с учетом числа ИЗА, дающих вклады в общую концентрацию в этой точке.

Перед расчетом из процесса нормирования исключаются ИЗА, дающие незначительный вклад в общую концентрацию в точке (меньше С_зн. )

Квота вклада источника определяется поэтапно и рассчитывается по следующему алгоритму:

Этап 1 – начальное значение квоты принимается равным C_зн

C_p =P_зн,

где С_зн =N\K_общ ;

K_общ – число вкладов в точке;

N – целевая концентрация в точке равная 1 ПДК; для зон санитарной охраны курортов, домов отдыха, зон городов и других территорий с повышенными требованиями к охране атмосферного воздуха N= 0.8 ПДК.

Этап 2 – определяется число нормируемых вкладов при квоте C_p , т.е. число вкладов больше расчетной квоты:

K_норм. =K_общ. -K_{ненорм.,}

где: K_норм. – число нормируемых вкладов (С_j > C_p );

K_{ненорм.} – число ненормируемых вкладов (С_j < C_p );

C_j – вклад j–го источника в суммарную концентрацию.

Этап 3 – определяется новое значение расчетной квоты (C’_p ):

C’_p =(N-C_{ненорм.} )/K_норм.

Где C_ненорм – сумма вкладов в долях ПДК, не превышающих текущего значения C_{p .}

Если число нормируемых источников K’_норм. по квоте C’_p меньше числа нормируемых вкладов K’_норм. по квоте C_p , то повторяем этапы 2 и 3, приняв C’_p =C_p .

В противном случае принимаем нормативную квоту C_н = C’_p .

Нормативное значение вклада ИЗА С_{i норм} в долях ПДК принимается равным нормативной квоте C_н :

C_{норм j .} =C_н

Мощность выброса ЗВ каждого существенного источника снижается пропорционально требуемому снижению вклада в точке:

Q_{норм j} =Q_j *(С_{норм j |} /С_j )

В качестве норматива мощности выброса X_j принимается наименьшее значение Q_{норм j} из рассчитанных во всех точках, где данный источник дает вклад в общую концентрацию.

Другой встречающийся в методической литературе метод расчетного определения ПДВ для группы источников носит название МРН-87 [42]. Суть его заключается в том, что все контрольные точки (где С>N) ранжируются в порядке убывания, после чего расчет начинается с точки 1. Всем источникам, определяющим заданный процент (95%-100%) загрязнения в этой точке устанавливается кратность снижения, равная С/N. На основе свойства линейности загрязнение в остальных точках пересчитывается. Если существуют точки, где превышение С>N сохранятся (там подключаются другие источники), то процедура повторяется. И так до тех пор, пока во всех точках не будет обеспечено соблюдение норматива N (ПДК).

Оба рассмотренных метода дают частное решение ранее рассмотренной системы линейных неравенств (С- существующая концентрация, N- ее нормативное значение, Х - ПДВ) для I контрольных точек (i=1,..,I). Отсутствие целевой функции не позволяет интерпретировать смысл рассчитанных таким образом ПДВ. Опыт работы с этими методами, реализованными в составе ПК ЭРА-ВОЗДУХ, показывает, что в различных случаях то один, то другой дает более выгодные для предприятия результаты в смысле максимального оставшегося после снижения суммарного выброса.

2 . Симплекс-метод

Симлекс-метод - это характерный пример итерационных вычислений. используемых при решении большинства оптимизационных задач.

В вычислительной схеме симплекс-метода реализуется упорядоченный процесс, при котором, начиная с некоторой исходной допустимой угловой точки (обычно начало координат), осуществляются последовательные переходы от одной допустимой экстремальной точки к другой до тех пор, пока не будет найдена точка, соответствующая оптимальному решению.

2.1 Общая характеристика симплекс метода

Симплекс метод - это универсальный метод для решения линейных систем уравнений или неравенств и линейного функционала [25].

Общая идея симплекс метода для решения ЗЛП (задачи линейного программирования) состоит в:

- умении находить начальный опорный план;

- наличии признака оптимальности опорного плана;

- умении переходить к нехудшему опорному плану.

Пусть ЗЛП представлена системой ограничений в каноническом виде:

.

Говорят, что ограничение ЗЛП имеет предпочтительный вид, если при неотрицательной правой части левая часть ограничений содержит переменную, входящую с коэффициентом, равным единице, а в остальные ограничения равенства - с коэффициентом, равным нулю.

Пусть система ограничений имеет вид

Сведем задачу к каноническому виду. Для этого прибавим к левым частям неравенств дополнительные переменные . Получим систему, эквивалентную исходной:

,

которая имеет предпочтительный вид

.

В целевую функцию дополнительные переменные вводятся с коэффициентами, равными нулю .

Пусть далее система ограничений имеет вид

Сведём её к эквивалентной вычитанием дополнительных переменных из левых частей неравенств системы. Получим систему

Однако теперь система ограничений не имеет предпочтительного вида, так как дополнительные переменные

входят в левую часть (при ) с коэффициентами, равными –1. Поэтому, вообще говоря, базисный план

не является допустимым. В этом случае вводится так называемый искусственный базис. К левым частям ограничений-равенств, не имеющих предпочтительного вида, добавляют искусственные переменные . В целевую функцию переменные , вводят с коэффициентом М в случае решения задачи на минимум и с коэффициентом -М для задачи на максимум, где М - большое положительное число. Полученная задача называется М-задачей, соответствующей исходной. Она всегда имеет предпочтительный вид.

Пусть исходная ЗЛП имеет вид

(2.1)

(2.2)

(2.3)

причём ни одно из ограничений не имеет предпочтительной переменной. М-задача запишется так:

(2.4)

(2.5)

, , (2.6)

Задача (2.4)-(2.6) имеет предпочтительный план. Её начальный опорный план имеет вид

Если некоторые из уравнений (2.2) имеют предпочтительный вид, то в них не следует вводить искусственные переменные.

Теорема. Если в оптимальном плане

(2.7)

М-задачи (2.4)-(2.6) все искусственные переменные , то план является оптимальным планом исходной задачи (2.1)-(2.3).

Для того чтобы решить задачу с ограничениями, не имеющими предпочтительного вида, вводят искусственный базис и решают расширенную М-задачу, которая имеет начальный опорный план

Решение исходной задачи симплексным методом путем введения искусственных переменных называется симплексным методом с искусственным базисом.

Если в результате применения симплексного метода к расширенной задаче получен оптимальный план, в котором все искусственные переменные , то его первые n компоненты дают оптимальный план исходной задачи.

Теорема. Если в оптимальном плане М-задачи хотя бы одна из искусственных переменных отлична от нуля, то исходная задача не имеет допустимых планов, т. е. ее условия несовместны.

Признаки оптимальности.

Теорема . Пусть исходная задача решается на максимум. Если для некоторого опорного плана все оценки неотрицательны, то такой план оптимален.

Теорема. Если исходная задача решается на минимум и для некоторого опорного плана все оценки являются неположительными, то такой план оптимален.

Для привидения системы ограничений неравенств к каноническому виду, необходимо в системе ограничений выделить единичный базис.

I. Ограничения вида «0»- ресурсные ограничения. Справа находится то что мы используем на производстве, слева - то что получаем. При таких ограничения вводят дополнительные переменные с коэффициентом «+1», образующие единичный базис. В целевую функцию эти переменные войдут с коэффициентом «0».

II. Ограничения вида «=». Часто бывает, что несмотря на то что ограничения имеют вид равенства, единичный базис не выделяется или трудно выделяется. В этом случае вводятся искусственные переменные для создания единичного базиса - Yi. В систему ограничений они входят с коэффициентом «1». а в целевую функцию с коэффициентом «M», стремящимся к бесконечности (при Fmin - «+M», при Fmax - «-M»).

III. Ограничения вида «0» - Плановые ограничения. Дополнительные переменные (X), несущие определенный экономический смысл - перерасход ресурсов или перевыполнение плана, перепроизводство, добавляются с коэффициентом «-1», в целевую функцию - с коэффициентом «0». А искусственные переменные (Y) как в предыдущем случае.

2.2 Алгоритм симплекс метода (первая симплекс таблица)

Пусть система приведена к каноническому виду.

X1+ q1,m+1 Xm+1 + …. + q1,m+n Xm+n = h1

X2+ q1,m+1 Xm+1 + …. + q1,m+n Xm+n = h1

X3+ q1,m+1 Xm+1 + …. + q1,m+n Xm+n = h1

……………………………………………………………….

Xm+ qm,m+1 Xm+1 + …. + qm,m+n Xm+n =hm

В ней m базисных переменных, k свободных переменных. m+k=n - всего переменных.

Fmin= C1X1+ C2X2+ C3X3+....+ CnXn

Все hi должны быть больше либо равны нулю, где i=1,2...m. На первом шаге в качестве допустимого решения принимаем все Xj=0 (j=m+1,m+2,...,m+k). При этом все базисные переменные Xi=Hi.

Для дальнейших рассуждений вычислений будем пользоваться первой симплекс таблицей (таблица1).

Таблица 1.

Первый столбец- коэффициенты в целевой функции при базисных переменных.

Второй столбец - базисные переменные.

Третий столбец - свободные члены (hi00).

Самая верхняя строка - коэффициенты при целевой функции.

Вторая верхняя строка - сами переменные, входящие в целевую функцию и в систему ограничений.

Основное поле симплекс метода - система коэффициентов из уравнения.

Последняя строка - служит для того, чтобы ответить на вопрос: «оптимален план или нет».

Индексная строка позволяет нам судить об оптимальности плана:

1. При отыскании Fmin в индексной строке должны быть отрицательные и нулевые оценки.

2. При отыскании Fmax в индексной строке должны быть нулевые и положительные оценки.

Переход ко второй итерации:

Для этого отыскиваем ключевой (главный) столбец и ключевую (главную) строку.

Ключевым столбцом является тот в котором находится наибольший положительный элемент индексной строки при отыскании Fmin или наименьший отрицательный элемент при отыскании Fmax.

Ключевой строкой называется та, в которой содержится наименьшее положительное частное от деления элементов столбца H на соответствующие элементы ключевого столбца.

На пересечении строки и столбца находится разрешающий элемент.

На этом этапе осуществляется к переходу к последующим итерациям.

Переход к итерациям:

1. Выводится базис ключевой строки, уступая место переменной из ключевого столбца со своим коэффициентом.

2. Заполняется строка вновь введенного базиса путем деления соответствующих элементов выделенной строки предыдущей итерации на разрешающий элемент.

3. Если в главной строке содержится нулевой элемент, то столбец, в котором находиться этот элемент переноситься в последующую итерацию без изменения.

4. Если в главном столбце имеется нулевой элемент, то строка, в которой он находиться переноситься без изменения в последующую итерацию.

5. Остальные элементы переносятся по формуле:

3 . Формализация поставленной задачи

Прежде всего, покажем, что характерные свойства МАД и известные сведения из теории линейных операторов позволяют экспертизу АОМ и установление ПДВ формализовать в виде двух зависимых математических задач.

Пусть суммарное загрязнение ВБ города отдельной примесью характеризуется функцией C(X,t) пространственных координат и времени. Загрязнение считается допустимым при C(X,t) < N(X) _, где N – норматив. Если C(X,t) > N(X), то необходимы АО мероприятия по достижению норматива. При их планировании из суммарного загрязнения атмосферы требуется выделить вклады Cj(X,t), j=1,…,J от J заданных источников, под которыми могут подразумеваться как отдельные трубы, аэрационные фонари и т.д., так и их совокупности, объединенные по различным признакам (по принадлежности к одному цеху, предприятию, ведомству по высоте выброса и т.д.). Пусть выбросы источников есть Q = (Q ₁ ,Q ₂ ,…,Q_J ). Предположим, что остальные параметры (высота, координаты и т.д.) в результате АОМ не изменяются. Тогда возникающую при экспертизе планов АОМ города задачу - определить изменение dC характеристик загрязнения ВБ по сравнению с базовым (предплановым) периодом – можно записать в виде:

dC = C₀ – C_P = A(Q) – A(X) (3.21)

где А- оператор модельной зависимости C от Q; величины C с индексом '₀ ' относятся к базовому периоду, а с индексом '_P ' - к ожидаемому после реализации запланированных АОМ. Заметим, что не только ожидаемый C_P , но и существующий уровень загрязнения C₀ , суммарное значение которого в некоторых точках промышленного города регулярно измеряется [78], требует в (1.21) модельного представления C₀ = A(Q₀ ), поскольку в общем случае методы контроля загрязнения не могут указать вклад конкретного источника в измеряемую величину. Соотношение (3.21) показывает, что формализация выделенной задачи сводится к построению и надлежащему применению оператора А, позволяющего переходить от выбросов к характеристикам загрязнения ВБ и различать заданный источник на фоне всех остальных.

Пользуясь линейностью модели ОНД-86 по выбросам источников можно представить загрязнение атмосферы в контрольных точках жилой зоны в виде линейной формы:

,

где - концентрация в i-ой точке, - выброс j-го источника, - вклад j-го источника в i-й точке, который в дальнейшем будем называть коэффициентом влияния. Отметим, что на практике (поскольку число источников и контрольных точек конечно) оператор А представляет собой матрицу, состоящую из коэффициентов влияния a_ij .

Санитарные требования приводят к системе линейных неравенств:

,

решение которой ищется на интервале

,

где - технологически обоснованный минимум выброса j-го источника.

На сегодняшний день в методической литературе описаны два метода нахождения решения поставленной задачи: МРН-87 [24] и метод равного квотирования [26]. Оба метода кратко рассмотрены в параграфе . Они дают частное решение поставленной системы неравенств из соображений удобства нахождения единственного решения. Однако любое предприятие заинтересовано в минимальных затратах, необходимых для установления нормативных выбросов . Для этой цели, к поставленной системе неравенств добавляем целевую функцию:

Где в общем случае - стоимость снижения на единицу выброса для j-го источника. В данном виде решение X_j дает минимум затрат на достижение нормативного загрязнения атмосферы.

В предположении

,

что эквивалентно .

Следовательно, задача сводится к поиску максимально допустимого по сумме сочетания выбросов X_j данного предприятия, позволяющего достичь нормативного загрязнения атмосферы. В связи с линейностью модели ОНД-86 по отношению к выбросам для поиска решения поставленной задачи может быть применён Симплекс-метод.

4. Программная реализация и пример практического применения

Для достижения поставленной задачи по разработке интегрированного в ПК ЭРА-ВОЗДУХ программного модуля расчета оптимальных ПДВ проведено изучение структуры файлов, в которых головной модуль передают данные расчетного блоку УПРЗА ЭРА. Разработаны соответствующие процедуры для автоматизированного чтения всех необходимых файлов.

4.1 Выбор загрязняющих веществ

После указания директории с данными для расчета, программа сканирует файлы в папке WORK по маске «htop*.ppp», таким образом выбирая загрязняющие вещества, для которых имеются начальные данные и возможно провести расчет.

Далее для каждого отмеченного вещества независимо от других будет производиться считывание значений и расчет ПДВ (Xj).

4.2 Обработка точек с повышенным загрязнением

procedure get_point (s:string;var countPoint:integer;var point_pdk:tExtArray);

Процедура извлекает значения N_i (ПДК), а так же количество точек I (i=1,..,I), в которых выбросы превышают N_i , из файла вида «htop*.ppp» в директории /WORK/.

Внутренняя структура htop*.ppp представляет собой текстовый файл, содержащий таблицу с данными о контрольных точках. Значение N_i содержатся в 8-ом столбце. Чаще всего N_i в точках равно 1 или 0.8 для особо охраняемых территорий (санатории, зоны отдыха).

4 .3 Обработка источников

procedure get_funnel(s:string; var countFunnel : integer ; var funnel_name : tsArray; var funnel_m:tExtArray;var funnel_min:tExtArray);

Процедура извлекает данные из файла вида «ist_*.txt», которые находятся в директории /DAT/. «ist_*.txt» - это текстовый файл, в котором в табличном виде представлена информация об источниках выбросов, в том числе:

· количество источников J;

· уникальный код источника;

· существующий выброс источника Qj;

· минимально возможный выброс Qjmin (не всегда указывается).

4.4 Обработка таблиц влияния источников на точки

procedure get_pointfunnel ( s : string; countPoint : integer; countfunnel : integer; funnel_name : tsArray; funnel_m : tExtArray; var pointfunnelx2 : tExtArrayx2; var point_cf : tExtArray);

Процедура извлекает коэффициенты влияния a_ij из файлов вида «10pd*.ppp», где «10pd*.ppp» - текстовый файл, содержащий отчет о результатах работы программы.

4.5 Применении симплекс-метода

procedure get_simplexsolve ( countPoint : integer; countFunnel : integer ; point_pdk : tExtArray; point_cf : tExtArray; funnel_m : tExtArray; funnel_min : tExtArray; pointfunnelx2 : tExtArrayx2; var x : tExtArray; var s_temp : string);

Процедура, используя данные расчетов программы «ЭРА-воздух», при помощи симплекс-метода рассчитывает оптимальные выбросы Xj для источников при заданных условиях.

4.6 Вывод полученных результатов

Результаты полученных вычислений выводятся на форму программы, а так же в файлы вида: «h_pd*.gpv», где * - это код вещества, для которого производился расчет.

4.7 Сравнение различных методов расчета ПДВ для реального предприятия

В расчете загрязнения атмосферы диоксидом азота на расчетном прямоугольнике 9 на 9 км с шагом 500 метров проведен расчет максимальных разовых концентраций, создаваемых 70-ю источниками выброса различного типа и высоты. В итоге получено поле максимальных концентраций, в которых есть области превышения норматива N, который в данном случае равен ПДК.

Картина загрязнения представлена на рисунке 1:

Рис 1. Загрязнение атмосферы в окрестности исследуемого предприятия.

Программный комплекс Эра-Воздух производит автоматический выбор точек превышения норматива N и позволяет найти расчетные значения ПДВ на основе методов МРН-87 и Метода равного квотирования.

При выполнении работы в программный комплекс добавлен новый модуль, в соответствии с алгоритмом описанным в пункте 5. Этот модуль позволяет рассчитать оптимальное значение ПДВ с использованием симплексного метода. Результаты расчетов ПДВ различными методами представлено ниже в таблицах 4.1, 4.2 и 4.3.

Таблица 4.1. Результаты расчета ПДВ (методом равного квотирования)

ПРИМЕСЬ=0301 Азот(IV) оксид (Азота диоксид)

Город :001 Кемерово, Объект : 0025

-----------------------------------------------------------------------------

| Код |Высота |Существую-|Минимально| Коэфф. | Расчетное |Кратность|

N | источника |источн.|щий выброс|возможный | норми- | значение |снижения |

п/п| выброса | м | г/с | выброс |рования | П Д В |выброса |

---|-----------|-------|----------|---г/с----|--------|----г/с----|---------|

15 00250010354 40.0 8.2390 0.0 0.485 3.9961 2.062

53 00250010944 40.0 8.2370 0.0 0.461 3.7941 2.171

Остальные источники не подлежат нормированию.

-----------------------------------------------------------------------------

В сумме по 0301 51.8545 0.0 0.832 43.1687 1.201

-----------------------------------------------------------------------------

Таблица 4.2. Результаты расчета ПДВ (метод МРН-87)

ПРИМЕСЬ=0301 Азот(IV) оксид (Азота диоксид)

Город :001 Кемерово, Объект :0025

-----------------------------------------------------------------------------

| Код |Высота |Существую-|Минимально| Коэфф. | Расчетное |Кратность|

N | источника |источн.|щий выброс|возможный | норми- | значение |снижения |

п/п| выброса | м | г/с | выброс |рования | П Д В |выброса |

---|-----------|-------|----------|---г/с----|--------|----г/с----|---------|

14 00250010353 33.4 1.0750 0.0 0.739 0.7942 1.354

15 00250010354 40.0 8.2390 0.0 0.739 6.0867 1.354

16 00250010356 60.0 2.1500 0.0 0.739 1.5884 1.354

24 00250010656 45.0 1.5000 0.0 0.739 1.1082 1.354

41 00250010892 60.0 9.0430 0.0 0.739 6.6807 1.354

52 00250010943 33.4 0.5280 0.0 0.739 0.3901 1.354

53 00250010944 40.0 8.2370 0.0 0.739 6.0852 1.354

54 00250010946 45.0 2.1470 0.0 0.739 1.5861 1.354

63 00250011198 24.4 0.5110 0.0 0.739 0.3775 1.354

66 00250011225 24.4 0.5890 0.0 0.739 0.4351 1.354

Остальные источники не подлежат нормированию.

----------------------------------------------------------------------------

В сумме по 0301 51.8545 0.0 0.829 42.9676 1.207

----------------------------------------------------------------------------

Таблица 4.3. Результаты расчета ПДВ (симплекс метод)

ПРИМЕСЬ=0301 Азот(IV) оксид (Азота диоксид)

Город :001 Кемерово, Объект :0025

---------------------------------------------------------

| Код |Существую-|Минимально| Расчетное | коэфф. |

| источника |щий выброс|возможный | значение | норми- |

| выброса | г/с | выброс | П Д В | рования |

|-----------|----------|---г/с----|----г/с----|---------|

|00250010353| 1.075000 | 0.000000 | 0.0196222 | 0.01825 |

|00250010943| 0.528000 | 0.000000 | 0.0000000 | 0.00000 |

|00250011225| 0.589000 | 0.000000 | 0.0000000 | 0.00000 |

- - - - - - - - - - - - - - - - - - - - - - - - - - - - -

| в сумме: 51.854470 49.682092 | 0.95811 |

---------------------------------------------------------

Решение задачи линейного программирования показывает, что максимально возможный выброс по заводу в целом можно обеспечить при достижении нормы загрязнения, если закрыть источники 1225 и 0943, а на источнике 0353 снизить выброс примерно в 50 раз. О реальности такого решения могут судить технологические службы предприятия.

Итоговая таблица сравнения трех методов:

Заключение

· Таким образом, в процессе выполнения дипломной работы проведен обзор существующей системы установления ПДВ для источников загрязнения атмосферы промышленных предприятий. При этом рассмотрена система обеспечения нормативных задач управления выбросами в атмосферу исходными данными. Исследована структура исходных данных их качество. Рассмотрены нормативные требования к загрязнению атмосферы населенных мест, которые накладывают ограничения на выбросы промышленных предприятий в атмосферу. Показано, что эти требования и свойство линейности нормативной модели расчета загрязнения атмосферы ОНД-86 по отношению к выбросам ИЗА позволяют представить процедуру нахождения расчетных значений ПДВ в виде задачи линейного программирования.

· Рассмотрены существующие в методической литературе способы расчета ПДВ в виде методов равного квотирования и МРН-87. Данные методы не имеют целевой функции и дают некоторые частные решения поставленной задачи.

· Изучен программный комплекс ЭРА-Воздух и форматы хранения и передачи данных между его модулями.

· Разработана процедура автоматического считывания исходных данных и результатов из ПК ЭРА-Воздух для полного обеспечения задачи линейного программирования исходными данными.

· Разработана программа расчета ПДВ на основе симплекс метода.

· Решена практическая задача по расчету ПДВ для одного из крупных предприятий г. Кемерово как встроенными в ПК ЭРА-ВОЗДУХ методами (МРН-87, равное квотирование), так и с использованием симплекс метода.

· Показано, что использование оптимизационного метода расчета ПДВ позволяет обеспечить нормативное загрязнение атмосферы при больших суммарных выбросах. В случае технологической приемлемости такого решения предприятие может существенно снизить платежи за сверхнормативный выброс в атмосферу.

Список литературы

1. Методическое пособие по расчету, нормированию и контролю выбросов загрязняющих веществ в атмосферный воздух. - СПб.: НИИ Атмосфера МПР РФ, 2002.

2. ОНД-86. Методика расчета концентраций в атмосферном воздухе вредных веществ, содержащихся в выбросах предприятий. Л.: Гидрометеоиздат, 1987.

3. Постановление Правительства Российской Федерации N182 от 2 марта 2000 г. «О порядке установления и пересмотра экологических и гигиенических нормативов качества атмосферного воздуха, предельно допустимых уровней физических воздействий на атмосферный воздух и государственной регистрации вредных (загрязняющих) веществ и потенциально опасных веществ». М., 2000.

4. Постановление Правительства Российской Федерации от 2 марта 2000 г. N 183 «О нормативах выбросов вредных (загрязняющих) веществ в атмосферный воздух и вредных физических воздействий на него». М., 2000.

5. Федеральный Закон «Об охране окружающей среды». М., 2002.

6. Рязанов В.А. О критериях и методах обоснования максимально допустимых концентраций атмосферных загрязнений в СССР.- В кн.: Предельно допустимые концентрации атмосферных загрязнений. Вып.8. - М.: Медицина, 1964, с. 5-21.

7. Перечень и коды веществ, загрязняющих атмосферный воздух. СПб., 2000.

8. Беспамятнов Г.П., Кротов Ю.А. Предельно допустимые концентрации химических веществ в окружающей среде - Л.: «Химия», 1985.

9. Пинигин М.А. Значение вероятностного подхода при решении вопросов гигиенического регламентирования атмосферных загрязнений. В кн. ”Медицинские проблемы охраны окружающей среды”. М.: 1981, с.95-102.

10. Берлянд М.Е. Прогноз и регулирование загрязнения атмосферы. –Л.: Гидрометеоиздат, 1985, 272с.

11. Безуглая Э.Ю. Мониторинг состояния загрязнения атмосферы в городах. Л.: Гидрометеоиздат, 1986, 200с.

12. Безуглая Э.Ю., Ковалевский А.Г., Расторгуева Г.П. Особенности распределения промышленных примесей в атмосфере городов различных типов. Тр. ГГО, вып. 467, 1983, с.81-87.

13. Перечень методик выполнения измерений концентраций загрязняющих веществ в выбросах промышленных предприятий СПб., 2001.

14. Перечень документов по расчету выделений (выбросов) загрязняющих веществ в атмосферный воздух, действующих в 2001-2002 годах. СПб., 2001.

15. Инструкция по инвентаризации выбросов загрязняющих веществ в атмосферу. Л., 1990.

16. СанПиН 2.1.6.1032-01 «Гигиенические требования к обеспечению качества атмосферного воздуха населенных мест». М., 2001.

17. .Атмосферная турбулентность и моделирование распространения примесей /под.ред. Ньистадта Ф.Т.М., Ван-Допа Х.- Л.: Гидрометеоиздат, 1985,-350 c.

18. Пененко В.В., Алоян А.Е. Модели и методы для задач охраны окружающей среды. -Новосибирск.: Наука, 1985.-256с.

19. Постановление Совета Министров РСФСР. Об утверждении на 1991 год нормативов за выбросы загрязняющих веществ в природную среду и порядка их применения./9 января 1991г. N 13 /. Собрание постановлений правительства РСФСР. -М.: N9, 1991.

20. Hanna S.R. Review of Atmospheric Diffusion Models for Regulatory Application.- WMO Tecnical Notes, No.177, 1982-42p.

21. Методы анализа загрязнений воздуха./Дугов Ю.С., Беликов А.Б., Дьяков Г.А., Тульчинский В.М.-М.: Химия, 1984,-384 с.

22. Вольберг Н.Ш., Егорова Е.Д., Кузьмина Т.А. Метрологические характеристики фотометрических методов анализа загрязнения атмосферы. - Тр. ГГО, 1982, No 450.c.107-111.

23. Израэль Ю.А, Гасилина Н.К., Ровинский Ф.Я. Система наблюдений и контроля загрязнения природной среды в СССР.- Метеорология и гидрология, 1978, No 10, c.5-12.

24. Методика расчета нормативов допустимых выбросов загрязняющих веществ в атмосферу для групп источников. МРН-87. - М., Госкомгидромет, Институт прикладной геофизики. 1987. -30с.

25. Рекомендации по определению допустимых вкладов в загрязнение атмосферы выбросов загрязняющих веществ предприятиями с использованием сводных расчетов загрязнения воздушного бассейна города (региона) выбросами промышленности и автотранспорта. СПб., 1999.-97с.

26. Васильев Ф.П. Методы решения экстремальных задач. М: Наука, 1980.-518с.

Вспомогательные указатели

Перечень сокращений

ЗВ - загрязняющее (вредное) вещество

ИЗА - источник загрязнения атмосферы

ПДВ - предельно допустимый выброс (допустимый выброс)

СЗЗ - санитарно-защитная зона

ПДК_р - максимальная разовая предельно допустимая концентрация загрязняющего вещества в атмосферном воздухе населенных мест

ПДК_с - среднесуточная предельно допустимая концентрация загрязняющего вещества в атмосферном воздухе населенных мест

ОБУВ - ориентировочный безопасный уровень воздействия загрязняющих веществ в атмосферном воздухе населенных мест

ГВС - газовоздушная смесь

ГОУ - газоочистная установка

ОНД - общесоюзный нормативный документ

НМУ - неблагоприятные метеорологические условия

УПРЗА - унифицированная программа расчета загрязнения атмосферы

Приложения

Unit1.pas

unit Unit1;

interface

uses

Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,

Dialogs, ShellAPI, ShlObj, StdCtrls, Buttons, CheckLst,Masks,inifiles,

ComCtrls,simplex, Menus;

const MyDecimalSeparator='.';

type

tsArray = array of string;

tExtArrayx2 = array of tExtArray;

TForm1 = class(TForm)

Edit1: TEdit;

GroupBox1: TGroupBox;

CheckListBox1: TCheckListBox;

Label1: TLabel;

BitBtn1: TBitBtn;

Button3: TButton;

Memo1: TMemo;

SpeedButton1: TSpeedButton;

CheckBox1: TCheckBox;

SpeedButton2: TSpeedButton;

SpeedButton3: TSpeedButton;

SpeedButton4: TSpeedButton;

procedure FormCreate(Sender: TObject);

procedure FormActivate(Sender: TObject);

procedure BitBtn1Click(Sender: TObject);

procedure Button3Click(Sender: TObject);

procedure N2Click(Sender: TObject);

procedure SpeedButton1Click(Sender: TObject);

procedure SpeedButton2Click(Sender: TObject);

procedure SpeedButton3Click(Sender: TObject);

procedure SpeedButton4Click(Sender: TObject);

private

{ Private declarations }

public

{ Public declarations }

end;

var

Form1: TForm1;

dir_path:string;

IniFile: TIniFile;

implementation

{$R *.dfm}

//запись в ini файл

procedure SaveIni(s:string);

var

IniPath: string;

FileName: string;

begin

GetDir(0,IniPath);

FileName:=IniPath+'\sav.ini';

IniFile:=TIniFile.Create(FileName);

Inifile.WriteString('patch','dir',s);

IniFile.Free;

end;

//чтение ini файла

function ReadIni:string;

var

IniPath: string;

FileName: string;

s:string;

begin

GetDir(0,IniPath);

FileName:=IniPath+'\sav.ini';

IniFile:=TIniFile.Create(FileName);

ReadIni:=Inifile.ReadString('patch','dir',s);

IniFile.Free;

end;

//--------- Удаляет пробел или запятую с краёв строки --------------------------

Function DelSpaceAndCap(s:string):string;

begin

while pos(copy(s,1,1),' ')<>0 do delete(s,1,1);

while pos(copy(s,length(s),1),' ')<>0 do delete(s,length(s),1);

result:=s;

end;

//--------- вырезает из строки имя ---------------------------------------------

Function ReturnSubString(Var s:string):string;

var

position,i : integer;

begin

s:=DelSpaceAndCap(s);

position:=0;

for i:=1 to length(s)-1 do

if (pos(copy(s,i,1),' ')<>0) and (position=0) then

position:=i;

if position=0 then begin

result:=s;

s:='';

end else begin

result := DelSpaceAndCap(copy(s,1,position));

Delete(s,1,position);

s:=DelSpaceAndCap(s);

end;

end;

//вывод ограничений

//==============================================================================

procedure vv(a:real;mas:tExtArray; Sign: TOperation);

var

i:integer;

s,s2,s3:string;

begin

s:=floattostr(mas[0]);

for i:=1 to length(mas)-1 do

s:=s+'+'+floattostr(mas[i]);

if Sign=less then s2:=' < ';

if Sign=Greater then s2:=' > ';

if Sign=Equal then s2:=' = ';

form1.memo1.lines.Add(s+s2+floattostr(a));

end;

//==============================================================================

//==============================================================================

//замена в строке всех вхождений одной подстроки на другую

function StrReplace(Str, Str1, Str2 : string):string;

var

p, L : integer;

s:string;

begin

s:=str;

L:=length(str1);

repeat

p:=pos(str1, s);

if p>0 then begin

Delete(s,p,L);

insert(str2, s, P);

end;

until P = 0;

StrReplace:=s;

end;

//==============================================================================

//==============================================================================

//========================= считывание таблиц влияния таблиц источников на точки

procedure get_pointfunnel(s:string;countPoint:integer;countfunnel:integer;funnel_name:tsArray;funnel_m:tExtArray;

var pointfunnelx2:tExtArrayx2; var point_cf:tExtArray);

var

h:textfile;

k,m:integer;

s_temp,s_temp2,s_temp3:string;

flag:boolean;

begin

SetLength(PointFunnelx2,countPoint,countFunnel);

SetLength(point_cf,countPoint);

for k:=0 to countPoint-1 do begin

point_cf[k]:=0;

for m:=0 to countFunnel-1 do

PointFunnelx2[k,m]:=0;

end;

AssignFile(h,dir_path+'\RESULT\'+'10pd'+s+'.ppp');

reset(h);

for k:=1 to 22 do readln(h,s_temp);

s_temp:=StrReplace(s_temp,'|',' ');

s_temp2:=s_temp;

for m:= 0 to CountPoint-1 do begin //общий цикл

flag:=true;

while flag do begin

if ReturnSubString(s_temp2)='Фоновая' then begin

point_cf[m]:=strtofloat(copy(s_temp,pos('%',s_temp)-4,4));

end else begin

s_temp3:=ReturnSubString(s_temp);

s_temp3:=ReturnSubString(s_temp);

s_temp3:=ReturnSubString(s_temp);

for k:=1 to 6 do s_temp2:=ReturnSubString(s_temp);

//showmessage(s_temp2);

for k:=0 to countFunnel-1 do

if s_temp3=copy(funnel_name[k],8,4) then

PointFunnelx2[m,k]:=strtofloat(s_temp2);//*funnel_m[k];

end;

readln(h,s_temp);

s_temp:=StrReplace(s_temp,'|',' ');

s_temp2:=s_temp;

if ReturnSubString(s_temp2)='В' then flag:=false;

end;

for k:=1 to 16 do readln(h,s_temp);

s_temp:=StrReplace(s_temp,'|',' ');

s_temp2:=s_temp;

end;

closefile(h);

end;

//==============================================================================

//==============================================================================

//========================================================= получение источников

procedure get_funnel(s:string; var countFunnel:integer;var funnel_name:tsArray;

var funnel_m:tExtArray;var funnel_min:tExtArray);

var

h,h2 : textfile;

index_funnel : integer;

i,j : integer;

s_temp,s_temp2:string;

begin

AssignFile(h,dir_path+'\DAT\'+'ist_'+s+'.txt');

reset(h);

index_funnel:=-11;

while s_temp<>'endI' do begin //чтение файла (установка размера массива)

readln(h,s_temp);

inc(index_funnel);

end;

closefile(h);

CountFunnel:=index_funnel;

setLength(funnel_m,CountFunnel);

setLength(funnel_min,CountFunnel);

setLength(funnel_name,CountFunnel);

for i:=0 to countFunnel-1 do begin

funnel_m[i]:=0;

funnel_min[i]:=0;

funnel_name[i]:='';

end;

AssignFile(h2,dir_path+'\DAT\'+'ist_'+s+'.txt');

reset(h2);

for j:=1 to 9 do

readln(h2,s_temp);

for i:= 0 to CountFunnel-1 do begin

readln(h2,s_temp);

funnel_name[i]:=ReturnSubString(s_temp);

for j:=1 to 14 do

s_temp2:=ReturnSubString(s_temp);

funnel_m[i]:=strtofloat(ReturnSubString(s_temp));

if DelSpaceAndCap(s_temp)<>'' then

funnel_min[i]:=strtofloat(DelSpaceAndCap(s_temp))

else funnel_min[i]:=0;

end;

closefile(h2);

end;

//==============================================================================

//==============================================================================

//============================================================= получение точек

procedure get_point (s:string;var countPoint:integer;var point_pdk:tExtArray);

var

index_point : integer;

i,j : integer;

h,h2 : textfile;

s_temp : string;

begin

index_point:=-2; // переменная для подсчета кол-ва точек

AssignFile(h,dir_path+'\WORK\'+'htop'+s+'.ppp');

reset(h);

while s_temp<>'000' do begin//чтение файла (установка размера массива)

readln(h,s_temp);

inc(index_point);

end;

closefile(h);

CountPoint:=index_point;

setLength(point_pdk,countPoint);

for i:=0 to countPoint-1 do

point_pdk[i]:=0; //зануление

AssignFile(h2,dir_path+'\WORK\'+'htop'+s+'.ppp');

reset(h2);

readln(h2,s_temp);

for i:= 0 to countPoint-1 do begin

readln(h2,s_temp);

for j:=1 to 8 do

point_pdk[i]:=strtofloat(ReturnSubString(s_temp));

end;

closefile(h2);

end;

//==============================================================================

//==============================================================================

//=========================================== решение при помощи симплекс метода

procedure get_simplexsolve(countPoint:integer;countFunnel:integer;point_pdk:tExtArray;

point_cf:tExtArray;funnel_m:tExtArray;funnel_min:tExtArray;

pointfunnelx2:tExtArrayx2;var x:tExtArray;var s_temp:string);

var

mas_temp : tExtArrayx2;

i,j : integer;

sim : TSimplex;

L : tExtArray;

begin

setLength(mas_temp,countFunnel,countFunnel);

setLength(L,countFunnel);

setLength(x,countFunnel);

for i:=0 to countFunnel-1 do

for j:=0 to countFunnel-1 do begin

if i=j then mas_temp[i,j]:=1 else mas_temp[i,j]:=0;

L[j]:=1;

end;

Sim:=TSimplex.Create(L,true);

for i:=0 to countPoint-1 do begin

//showmessage(vv(point_pdk[i],pointfunnelx2[i]));

Sim.AddCons(point_pdk[i],pointfunnelx2[i],less);

if form1.CheckBox1.Checked then vv(point_pdk[i],pointfunnelx2[i],less);

end;

for i:=0 to countFunnel-1 do begin

Sim.AddCons(funnel_m[i],mas_temp[i],less);

if funnel_min[i]>0 then begin

Sim.AddCons(funnel_min[i],mas_temp[i],Greater);

if form1.CheckBox1.Checked then vv(funnel_min[i],mas_temp[i],Greater);

end;

end;

if (Sim.Solve=SIMPLEX_DONE) then begin

s_temp:='решение найдено';

x:=Sim.GetSolution;

end

else s_temp:='Решения не существует';

end;

//==============================================================================

//==============================================================================

//==================================================== общий модуль для подсчета

procedure TForm1.Button3Click(Sender: TObject);

var

s,s_temp,ss : string;

countPoint : integer;

countfunnel : integer;

point_pdk : tExtArray;

point_cf : tExtArray;

funnel_m : tExtArray;

funnel_min : tExtArray;

funnel_name : tsArray;

pointfunnelx2 : tExtArrayx2;

i,j : integer;

x : tExtArray;

empty : boolean;

h : textfile;

funnelSumM,sumX:real;

begin

funnelSumM:=0;

sumX:=0;

memo1.Clear;

for i:=0 to checkListBox1.Items.Count-1 do begin

if CheckListBox1.Checked[i] then begin

application.ProcessMessages;

s:=checklistbox1.Items.Strings[i];

s:=returnSubString(s);

application.ProcessMessages;

get_point (s,countPoint,point_pdk);

get_funnel(s,countFunnel,funnel_name,funnel_m,funnel_min);

get_pointfunnel(s,countPoint,countfunnel,funnel_name,funnel_m,pointfunnelx2,point_cf);

get_simplexsolve(countPoint,CountFunnel,point_pdk,point_cf,funnel_m,funnel_min,pointfunnelx2,x,s_temp);

AssignFile(h,dir_path+'\RESULT\'+'h_pd'+s+'.gpv');

rewrite(h);

if s_temp='решение найдено' then begin

memo1.lines.Add('');

memo1.lines.Add(' Результаты расчета ПДВ (симплекс метод):');

memo1.lines.Add(' ПРИМЕСЬ='+s);

memo1.lines.Add('');

memo1.lines.Add('---------------------------------------------------------');

memo1.lines.Add('| Код |Существую-|Минимально| Расчетное | коэфф. |');

memo1.lines.Add('| источника |щий выброс|возможный | значение | норми- |');

memo1.lines.Add('| выброса | г/с | выброс | П Д В | рования |');

memo1.lines.Add('|-----------|----------|---г/с----|----г/с----|---------|');

writeln(h,'');

writeln(h,' Результаты расчета ПДВ (симплекс метод):');

writeln(h,' ПРИМЕСЬ='+s);

writeln(h,'');

writeln(h,'---------------------------------------------------------');

writeln(h,'| Код |Существую-|Минимально| Расчетное | коэфф. |');

writeln(h,'| источника |щий выброс|возможный | значение | норми- |');

writeln(h,'| выброса | г/с | выброс | П Д В | рования |');

writeln(h,'|-----------|----------|---г/с----|----г/с----|---------|');

empty:=true;

for j:=0 to countFunnel-1 do begin

funnelSumM:=FunnelSumM+funnel_m[j];

sumX:=SumX+x[j];

if abs(x[j]-funnel_m[j])>0.0000001 then

begin

ss:='|'+funnel_name[j]+'| '+FloatToStrF(funnel_m[j],ffFixed,1000,6)+' | '+FloatToStrF(funnel_min[j],ffFixed,1000,6);

ss:=ss+' | '+FloatToStrF(x[j],ffFixed,1000,7)+' | '+FloatToStrF(x[j]/funnel_m[j],ffFixed,1000,5)+' |';

memo1.lines.Add(ss);

writeln(h,ss);

empty:=false;

end;

end;

ss:='| в сумме: '+FloatToStrF(funnelSumM,ffFixed,1000,6)+' ';

ss:=ss+FloatToStrF(sumX,ffFixed,1000,6)+' | '+ FloatToStrF(sumX/funnelSumM,ffFixed,1000,5)+' |';

if empty then begin

memo1.lines.Add('| Нет выбросов для снижения |');

writeln(h,'| Нет выбросов для снижения |');

end;

if not empty then begin

memo1.lines.Add('- - - - - - - - - - - - - - - - - - - - - - - - - - - - -');

memo1.lines.Add(ss);

writeln(h,'- - - - - - - - - - - - - - - - - - - - - - - - - - - - -');

writeln(h,ss);

end;

memo1.lines.Add('---------------------------------------------------------');

memo1.lines.Add('');

memo1.lines.Add('');

writeln(h,'---------------------------------------------------------');

writeln(h,'');

writeln(h,'');

end else begin

memo1.lines.Add('');

memo1.lines.Add(' Результаты расчета ПДВ (симплекс метод):');

memo1.lines.Add(' ПРИМЕСЬ='+s);

memo1.lines.Add('');

memo1.lines.Add('---------------------------------------------------------');

memo1.lines.Add('| Решение не найдено |');

memo1.lines.Add('---------------------------------------------------------');

writeln(h,'');

writeln(h,' Результаты расчета ПДВ (симплекс метод):');

writeln(h,' ПРИМЕСЬ='+s);

writeln(h,'');

writeln(h,'---------------------------------------------------------');

writeln(h,'| Решение не найдено |');

writeln(h,'---------------------------------------------------------');

end;

closefile(h);

end;

// closefile(h);

end;

end;

//==============================================================================

//поиск файла по маске

procedure FindFiles(StartFolder, Mask: string; List: TStrings;

ScanSubFolders: Boolean = True);

var

SearchRec: TSearchRec;

FindResult: Integer;

begin

List.BeginUpdate;

try

StartFolder := IncludeTrailingBackslash(StartFolder);

FindResult := FindFirst(StartFolder + '*.*', faAnyFile, SearchRec);

try

while FindResult = 0 do

with SearchRec do

begin

if (Attr and faDirectory) <> 0 then

begin

if ScanSubFolders and (Name <> '.') and (Name <> '..') then

FindFiles(StartFolder + Name, Mask, List, ScanSubFolders);

end

else

begin

if MatchesMask(Name, Mask) then begin

List.Add(copy(Name,5,4));

//showmessage(StartFolder + Name);

end;

end;

FindResult := FindNext(SearchRec);

end;

finally

FindClose(SearchRec);

end;

finally

List.EndUpdate;

end;

end;

procedure TForm1.FormCreate(Sender: TObject);

begin

DecimalSeparator:=MyDecimalSeparator;

end;

procedure TForm1.FormActivate(Sender: TObject);

begin

dir_path:=ReadIni;

edit1.Text:=dir_path;

{--}

end;

procedure TForm1.BitBtn1Click(Sender: TObject);

var

h,h2:textfile;

i,j,k,n:integer;

s_temp:string;

s: array of array of string;

begin

dir_path:=edit1.Text;

checklistbox1.Items.Clear;

i:=0;

AssignFile(h,dir_path+'\WORK\activ2.txt');

reset(h);

//readln(h,s_temp);

while not EOF(h) do begin//чтение файла (установка размера массива)

readln(h,s_temp);

inc(i);

end;

closefile(h);

setlength(s,i,2);

AssignFile(h2,dir_path+'\WORK\activ2.txt');

reset(h2);

for j:=0 to i-1 do begin

readln(h2,s_temp);

s[j,0]:=copy(s_temp,24,4);

s[j,1]:=copy(s_temp,30,55);

end;

closefile(h2);

FindFiles(dir_path, 'htop*.ppp', checklistbox1.items, true);

n:=checklistbox1.items.Count-1;

for j:=0 to n do begin

for k:=0 to i-1 do begin

//showmessage(s[k,0]+' -| ');

if checklistbox1.items[0]=s[k,0] then begin

//showmessage(s[j,0]+' | '+s[j,1]);

checklistbox1.items.Delete(0);

checklistbox1.items.Add(s[k,0]+' '+s[k,1]);

end;

end;

end;

end;

procedure TForm1.N2Click(Sender: TObject);

var

TitleName : string;

lpItemID : PItemIDList;

BrowseInfo : TBrowseInfo;

DisplayName : array[0..MAX_PATH] of char;

TempPath : array[0..MAX_PATH] of char;

begin

FillChar(BrowseInfo, siCeof(TBrowseInfo), #0);

BrowseInfo.hwndOwner := Form1.Handle;

BrowseInfo.psCDisplayName := @DisplayName;

TitleName := 'Please specify a directory';

BrowseInfo.lpsCTitle := PChar(TitleName);

BrowseInfo.ulFlags := BIF_RETURNONLYFSDIRS;

lpItemID := SHBrowseForFolder(BrowseInfo);

if lpItemId <> nil then

begin

SHGetPathFromIDList(lpItemID, TempPath);

edit1.Text:=TempPath;

GlobalFreePtr(lpItemID);

end;

//showmessage(tempPath);

dir_path:=tempPath;

//FindFiles(tempPath, 'htop*.ppp', checkmemo1.lines, true); //старая версия

SaveIni(dir_path);

end;

procedure TForm1.SpeedButton1Click(Sender: TObject);

var

TitleName : string;

lpItemID : PItemIDList;

BrowseInfo : TBrowseInfo;

DisplayName : array[0..MAX_PATH] of char;

TempPath : array[0..MAX_PATH] of char;

begin

FillChar(BrowseInfo, siCeof(TBrowseInfo), #0);

BrowseInfo.hwndOwner := Form1.Handle;

BrowseInfo.psCDisplayName := @DisplayName;

TitleName := 'Please specify a directory';

BrowseInfo.lpsCTitle := PChar(TitleName);

BrowseInfo.ulFlags := BIF_RETURNONLYFSDIRS;

lpItemID := SHBrowseForFolder(BrowseInfo);

if lpItemId <> nil then

begin

SHGetPathFromIDList(lpItemID, TempPath);

edit1.Text:=TempPath;

GlobalFreePtr(lpItemID);

end;

//showmessage(tempPath);

dir_path:=tempPath;

//FindFiles(tempPath, 'htop*.ppp', checkmemo1.lines, true); //старая версия

SaveIni(dir_path);

end;

procedure TForm1.SpeedButton2Click(Sender: TObject);

var

i:integer;

begin

for i:=0 to checklistbox1.Items.Count-1 do

checklistbox1.Checked[i]:=true;

end;

procedure TForm1.SpeedButton3Click(Sender: TObject);

var

i:integer;

begin

for i:=0 to checklistbox1.Items.Count-1 do

checklistbox1.Checked[i]:=false;

end;

procedure TForm1.SpeedButton4Click(Sender: TObject);

var

i:integer;

begin

for i:=0 to checklistbox1.Items.Count-1 do

if checklistbox1.Checked[i] then checklistbox1.Checked[i]:=false

else checklistbox1.Checked[i]:=true;

end;

end.

Simplex.pas

unit simplex;

interface

const

SIMPLEX_DONE = 0; // оптимизация успешно завершена

SIMPLEX_NO_SOLUTION = 1; // задача не имеет решения (не удается найти базис)

SIMPLEX_NO_BOTTOM = 2; // решения нет, т.к. линейная форма не ограничена снизу

SIMPLEX_NEXT_STEP = 3; // для получения решения нужно сделать еще хотя бы один шаг

MAX_VAL = 0.1e-12; //точность (значение, удовлетворяющее -MAX_VAL < X < MAX_VAL считается нулем)

type

TOperation = (Equal,Less,Greater);

TExtArray = array of extended;

TConstrain = record

A : TExtArray;

B : extended;

Sign : TOperation;

isT : boolean;

end;

TSimplex = class

M,N : integer; { M - число строк, N - число столбцов}

RealN : integer; {реальное число переменных, изначально вошедших в задачу}

Cons : array of TConstrain;

C : TExtArray;

L : extended;

Basis : array of integer;

Max : boolean; { направление оптимизации: минимизация или максимизация }

Constructor Create(_C:TExtArray; MaximiCe:boolean=false);

Constructor CreateBasis(const Simplex:TSimplex);

Constructor Copy(const Simplex:TSimplex);

Procedure AddCons(_B:extended; _A:TExtArray; Sign:TOperation);

Procedure SetAllLengths(Len:integer);

Function SimplexStep:integer;

Function CheckBasis:boolean;

Function FoundInBasis(num:integer): integer;

Function DoPrec(num:extended): extended;

Procedure NormaliCe;

Procedure MulString(Number:integer; Value:extended);

Procedure AddString(Num1,Num2:integer; Value:extended); {суммирование строки 1 со строкой 2, домноженной на коэффициент Value }

Function Solve:integer;

Function GetMin:extended;

Function GetSolution:TExtArray;

Destructor Free;

end;

TIntSimplex = class(TSimplex)

// CurX : TExtArray;

//CurL : extended;

// CurFound : boolean;

Constructor Create(_C:TExtArray; MaximiCe:boolean=false);

// Procedure DelLastCons;

Function IntSolve:integer;

Function GetIntMin:extended;

Function IsInteger(value:extended):boolean;

Function GetIntSolution:TExtArray;

// Function SearchCons(_B:extended;_A:TExtArray):integer;

end;

implementation

uses Math;

{ TSimplex }

Function TSimplex.DoPrec(num:extended): extended;

begin

if ((num < MAX_VAL) and (num > -MAX_VAL)) then

num := 0;

Result := num;

end;

procedure TSimplex.AddCons(_B: extended; _A: TExtArray; Sign: TOperation);

var

j : integer;

begin

if (Length(_A)>N) then SetAllLengths(Length(_A));

inc(M);

SetLength(Cons,M);

//if ((_B=0) and (Sign=Less)) then Sign:=Equal; //???

Cons[M-1].B:=_B;

Cons[M-1].Sign:=Sign;

SetLength(Cons[M-1].A,N);

for j:=0 to Length(_A)-1 do Cons[M-1].A[j]:=_A[j];

if Length(_A)<N then for j:=Length(_A) to N-1 do Cons[M-1].A[j]:=0;

end;

{суммирование строки 1 со строкой 2, домноженной на коэффициент Value }

procedure TSimplex.AddString(Num1, Num2: integer; Value: extended);

var

j : integer;

begin

for j:=0 to N-1 do Cons[Num1].A[j]:=Cons[Num1].A[j]+Cons[Num2].A[j]*Value;

Cons[Num1].B:=Cons[Num1].B+Cons[Num2].B*Value;

end;

function TSimplex.CheckBasis: boolean;

var

i,j,k : integer;

f : boolean;

begin

SetLength(Basis,M);

for i:=0 to M-1 do Basis[i]:=-1;

for j:=0 to N-1 do begin

f:=true;

k:=-1;

i:=0;

while (f and (i<M)) do begin

if ((Cons[i].A[j]<>0) and (Cons[i].A[j]<>1)) then f:=false;

if (Cons[i].A[j]=1) then begin

if (k=-1) then k:=i

else f:=false;

end;

inc(i);

end;

if (f and (k<>-1)) then Basis[k]:=j;

end;

f:=true;

for i:=0 to M-1 do f:=f and (Basis[i]<>-1);

Result:=f;

end;

constructor TSimplex.Create(_C: TExtArray; MaximiCe:boolean);

var

j : integer;

begin

N:=Length(_C);

RealN := N;

M:=0;

SetLength(C,N);

Max:=MaximiCe;

if (not MaximiCe) then for j:=0 to N-1 do C[j]:=-_C[j]

else for j:=0 to N-1 do C[j]:=_C[j];

Max:=MaximiCe;

L := 0;

end;

constructor TSimplex.Copy(const Simplex: TSimplex);

var

i,j : integer;

begin

M:=Simplex.M;

N:=Simplex.N;

RealN := Simplex.RealN;

SetLength(Cons,M);

SetLength(Basis,M);

SetLength(C,N);

Max:=Simplex.Max;

for i:=0 to M-1 do begin

SetLength(Cons[i].A,N);

Basis[i]:=-1;

for j:=0 to N-1 do Cons[i].A[j]:=Simplex.Cons[i].A[j];

Cons[i].B:=Simplex.Cons[i].B;

Cons[i].Sign:=Simplex.Cons[i].Sign;

end;

for i:=0 to Simplex.N-1 do C[i]:=Simplex.C[i];

L := Simplex.L;

end;

constructor TSimplex.CreateBasis(const Simplex: TSimplex);

var

i,j : integer;

begin

M:=Simplex.M;

N:=Simplex.N;

RealN := Simplex.RealN;

L := 0;

SetLength(Cons,M);

SetLength(Basis,M);

SetLength(C,N);

for i:=0 to N-1 do C[i]:=0;

for i:=0 to M-1 do begin

SetLength(Cons[i].A,N);

for j:=0 to N-1 do Cons[i].A[j]:=Simplex.Cons[i].A[j];

Cons[i].B:=Simplex.Cons[i].B;

Cons[i].Sign:=equal;

Cons[i].isT := false;

end;

for i:=0 to M-1 do begin

if (Simplex.Basis[i]<>-1) then Basis[i]:=Simplex.Basis[i]

else begin

SetAllLengths(N+1);

for j:=0 to M-1 do Cons[j].A[N-1]:=0;

Cons[i].A[N-1]:=1;

Cons[i].isT := true;

C[N-1] := 0;

for j:=0 to Simplex.N-1 do C[j] := C[j] + Simplex.Cons[i].A[j];

L := L + Cons[i].B;

end;

end;

end;

destructor TSimplex.Free;

begin

SetLength(C,0);

SetLength(Basis,0);

SetLength(Cons,0);

M:=0;

N:=0;

RealN := 0;

end;

function TSimplex.GetMin: extended;

var

i : integer;

begin

if (Max) then

Result := -L

else

Result := L;

end;

function TSimplex.GetSolution: TExtArray;

var

Solution : TExtArray;

i,j : integer;

begin

SetLength(Solution,RealN);

for j:=0 to RealN-1 do begin

Solution[j]:=0;

i:=0;

while ((i<M) and (Basis[i]<>j)) do inc(i);

if ((Basis[i]=j) and (i<M)) then Solution[j]:=Cons[i].B;

end;

Result:=Solution;

end;

procedure TSimplex.MulString(Number: integer; Value: extended);

var

j : integer;

begin

for j:=0 to N-1 do Cons[Number].A[j]:=Cons[Number].A[j]*Value;

Cons[Number].B:=Cons[Number].B*Value;

end;

procedure TSimplex.NormaliCe;

var

i : integer;

begin

for i:=0 to M-1 do if (Cons[i].Sign<>Equal) then begin

SetAllLengths(N+1);

if (Cons[i].Sign=Greater) then Cons[i].A[N-1]:=-1

else Cons[i].A[N-1]:=1;

Cons[i].Sign := Equal;

end;

end;

procedure TSimplex.SetAllLengths(Len: integer);

var

i, j : integer;

OldN : integer;

begin

OldN:=N;

N:=Len;

SetLength(C,N);

for i:=0 to M-1 do SetLength(Cons[i].A,N);

if (OldN<N) then begin

for j:=OldN to N-1 do begin

C[j]:=0;

for i:=0 to M-1 do Cons[i].A[j]:=0;

end;

end;

end;

function TSimplex.FoundInBasis(num:integer): integer;

var

i:integer;

f:boolean;

begin

f := false;

i := 0 ;

while (not f and (i<M)) do

begin

f := (Basis[i] = num);

inc(i);

end;

if (f) then

Result := i-1

else

Result := -1;

end;

function TSimplex.SimplexStep: integer;

var

i,j : integer;

f,opt : boolean;

x,y : integer; //координаты опорного элемента

CurMax : extended;

temp : array of TConstrain;

tempC : TExtArray;

begin

opt := true;

CurMax := -1;

for i := 0 to N-1 do

begin

//проверка на разрешимость

if (C[i] > 0) then

begin

opt := false; //а это попутная проверка на оптимальность

if (C[i] > CurMax) then //а это поиск ведущего столбца (максимальный элемент в C[i])

begin

CurMax := C[i];

x := i;

end;

f := true;

for j := 0 to M-1 do

f := f and (Cons[j].A[i] < 0);

if (f) then

begin

Result := SIMPLEX_NO_BOTTOM;

exit;

end;

end;

end;

if (opt) then

Result := SIMPLEX_DONE

else

begin

//зная номер ведущего столбца, ищем номер ведущей строки

CurMax := MaxExtended; //на самом деле тут будем искать минимум, а не Max

for j := 0 to M-1 do

if (Cons[j].A[x] > 0) then //идем только по положительным элементам

if (Cons[j].B/Cons[j].A[x] < CurMax) then

begin

CurMax := Cons[j].B/Cons[j].A[x];

y := j;

end

else if (DoPrec(Cons[j].B/Cons[j].A[x] - CurMax) = 0) then

if (Cons[j].isT) then

y := j;

//сохраняем текущие значения

SetLength(temp, M);

for j := 0 to M-1 do

begin

SetLength(temp[j].A, N);

for i := 0 to N-1 do

temp[j].A[i] := Cons[j].A[i];

temp[j].B := Cons[j].B;

end;

SetLength(tempC, N);

for i := 0 to N-1 do

tempC[i] := C[i];

//делаем пересчет таблицы

//строка делиться на ведущий элемент

MulString(y, 1/Cons[y].A[x]);

//преобразование остальных элементов

for j := 0 to M-1 do

begin

if (j <> y) then

begin

for i := 0 to N-1 do

begin

Cons[j].A[i] := DoPrec(temp[j].A[i] - temp[j].A[x]*temp[y].A[i]/temp[y].A[x]);

end;

Cons[j].B := DoPrec(temp[j].B - temp[j].A[x]*temp[y].B/temp[y].A[x]);

end

else

begin

for i := 0 to N-1 do

Cons[j].A[i] := DoPrec(Cons[j].A[i]);

end;

end;

//и строка с коэффициентами функции

for i := 0 to N-1 do

begin

C[i] := DoPrec(tempC[i] - tempC[x]*temp[y].A[i]/temp[y].A[x]);

end;

Basis[y] := x;

//и сама функция:

L := DoPrec(L - tempC[x]*temp[y].B/temp[y].A[x]);

for i:= 0 to M-1 do

SetLength(temp[i].A, 0);

SetLength(temp, 0);

SetLength(tempC, 0);

Result := SIMPLEX_NEXT_STEP;

end;

end;

function TSimplex.Solve: integer;

var

i,j : integer;

Simplex : TSimplex;

f : boolean;

Step : integer;

cc : extended;

begin

//oldN := N;

NormaliCe;

f:=false;

if (not CheckBasis) then begin

Simplex:=TSimplex.CreateBasis(self);

Simplex.Solve;

f:=Simplex.GetMin<>0;

if (not f) then for i:=0 to M-1 do begin

for j:=0 to N-1 do Cons[i].A[j]:=Simplex.Cons[i].A[j];

Cons[i].B:=Simplex.Cons[i].B;

Cons[i].isT := false;

Basis[i]:=Simplex.Basis[i];

cc := C[Basis[i]];

for j:=0 to N-1 do

C[j] := DoPrec(C[j] - cc*Cons[i].A[j]);

L := DoPrec(L - cc*Cons[i].B);

end;

Simplex.Free;

end;

if (f) then Step:=SIMPLEX_NO_SOLUTION

else repeat

Step:=SimplexStep;

until (Step<>SIMPLEX_NEXT_STEP);

//SetAllLengths(OldN);

Result:=Step;

end;

{ TIntSimplex }

constructor TIntSimplex.Create(_C:TExtArray; MaximiCe:boolean=false);

begin

//CurFound:=false;

inherited;

end;

function TIntSimplex.GetIntMin: extended;

begin

Result:=GetMin;

end;

function TIntSimplex.GetIntSolution: TExtArray;

begin

Result:=GetSolution;

end;

function TIntSimplex.IsInteger(Value:extended):boolean;

begin

Result:=((Value=floor(Value)) or (Value=ceil(Value)));

end;

function TIntSimplex.IntSolve: integer;

var

i : integer;

OldN : integer;

FractCol : integer;

FractRow : integer;

TmpX : TExtArray;

TmpCons : TExtArray;

NewValue : extended;

begin

if (Solve=SIMPLEX_DONE) then begin

//if (not CurFound or ((Simplex.GetMin<CurL) and not Max) or ((Simplex.GetMin>CurL) and Max)) then begin

TmpX:=GetSolution;

i:=0;

while ((i<RealN) and IsInteger(TmpX[i])) do inc(i);

FractCol:=i;

if (FractCol<>RealN) then begin // если найдена хотя бы одна нецелая переменная

OldN:=N;

SetLength(TmpCons,N);

FractRow := FoundInBasis(FractCol);

for i := 0 to N-1 do

if (FoundInBasis(i) = -1) then

TmpCons[i] := Cons[FractRow].A[i] - Floor(Cons[FractRow].A[i])

else

TmpCons[i] := 0;

NewValue := Cons[FractRow].B - Floor(Cons[FractRow].B);

//if (Max) then

AddCons(NewValue, TmpCons, Greater);

//else

// AddCons(NewValue, TmpCons, Less);

Result := IntSolve;

SetAllLengths(OldN); // удаляем пустые столбцы в конце, если они есть

end

else begin // если полученное решение - целочисленное\

Result := SIMPLEX_DONE;

end;

//end;

end

else

Result:=SIMPLEX_NO_SOLUTION;

end;

end.