Бесплатный автореферат и диссертация по геологии на тему

Физико-химические основы моделирования гидротермальных систем

ВАК РФ 04.00.02, Геохимия

Автореферат диссертации по теме "Физико-химические основы моделирования гидротермальных систем"

р V Б ОА

1 О

На правах рукописи

Акинфиев Николай Николаевич

Физико-химические основы моделирования гидротермальных систем

04.00.02 - геохимия Автореферат

диссертации на соискание учёной степени доктора химических наук

Москва - 1995

Работа выполнена в Московской государственной геологоразведочной академии им. С. Орджоникидзе

Научный консультант: доктор хим. наук Рыженко Б.Н.

Официальные оппоненты: доктор хим. наук, профессор Кесслер Ю.М.

доктор хим. наук, профессор Малинин С. Д. доктор техн. наук, профессор Путилов A.B. Ведущая организация: Институт экспериментальной минералогии РАН

Защита состоится ¿С апреля 1995 г. в на заседании

диссертационного совета Д002.59.02 Института геохимии и аналитической химии им. В. И. Вернадского РАН по адресу Москва В-334, ул. Косыгина 19.

С диссертацией можно ознакомиться в библиотеке ГЕОХИ РАН. Автореферат разослан

Ученый секретарь диссертационного совета

Жидикова А.П.

г

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

Во-первых, это разработка алгоритмов и создание программ расчета равновесий в многокомпонентных гетерогенных системах. К моменту начала работы над диссертацией наиболее эффективные и универсальные алгоритмы, основанные на методе прямой минимизации функции Гиббса системы [Карпов И.Н.. 1976; Шваров Ю.В.. 1978]. были реализованы лишь на "больших" ЭВМ. При этом весьма актуальной была задача создания универсального алгоритма расчёта для персональных микро-ЭВМ с ограниченным объёмом оперативной памяти.

Во-вторых, это решение фундаментальных задач предсказания термодинамических свойств водных частиц при повышенных параметрах состояния. В настоящее время для этих целей широко применяется уравнение состояния Хелгесона-Киркхэма-Флауэрса (НКП [Не^еэоп е1 а1, 1981, 1988. 1992]. Использование этого уравнения позволяет с достаточной точностью описывать термодинамические свойства водных частиц в широком диапазоне внешних условий (до 1000°С и до 5 кбар). Однако индивидуальные параметры уравнения НЮ7 носят в значительной степени эмпирический характер, что значительно снижает предсказательный потенциал модели. Аналогичная проблема возникает при попытке описания комплексообразования в водных растворах при повышенных температурах. Для описания ионной ассоциации в этих условиях обычно используются различные варианты электростатической теории (Риоэз И.М., 1958; Ры-женко Б.Н.. 1981) также не свободные от эмпирических параметров. В этой связи весьма актуальной становится задача отыскания физически обоснованнго и по возможности лишённого подгоночных параметров подхода для описания поведения водных частиц в широком диапазоне температур и давлений.

И наконец, в-третьих, это получение надежной термодинамической информации для химических компонентов конкретных систем. Для исследуемых в настоящей работе Аэ11 '-Б11 —Н20 и БЬ1 '-Б1 !-С1-Н20 в литературе имелся обширный материал взаимно несогласованных экспериментальных данных по растворимости Аэ- и БЬ-содержащих фаз. Таким образом, актуальной являлась задача выбора наиболее надежных и внутренне согласованных экспериментальных данных, их новой интерпретации, и путем решения обратных термодинамических задач получения согласованных значений термодинамических параметров компонентов системы.

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

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

Разработаны методы решения обратных термодинамических задач, когда по имеющимся экспериментальным данным восстанавливаются согласованные значения термодинамических параметров компонентов системы, позволяющие затем проводить интерполяцию физико-химических свойств в широком диапазоне давлений и температур. Кроме получения новой термодинамической информации по компонентам систем Аэ11 ^Б11 -Н20 и БЬ111 -Б1 '-СИ-НгО научная новизна развиваемого подхода заключалась в последовательном использовании уравнения состояния ИКР для частиц водного раствора.

Предложена модель кинетики растворения некоторых минералов, позволяющая с новых позиций объяснить основные экспериментальные закономерности скорости их растворения.

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

Используемый подход впервые позволил перейти от рассмотрения отдельных реакций к расчету равновесия системы в целом. Расчёт растворимости галенита в условиях, соответствующих эволюции модельного H20-C02-fJaCl-H2S флюида показал, что изменение диэлектрической проницаемости при вскипании можно рассматривать как новый мощный фактор минералообразования.

Практическая значимость. Предложенные в работе подходы необходимы для понимания процессов, сопровождающих мобилизацию, транспорт и отложение рудного вещества и, следовательно, построения количественных моделей рудоотложения. Широкое практическое использование результатов предполагается при исследовании гидротермальных процессов (геохимия гидротермального минералообразования, высокотемпературная химическая технология), для организации баз данных и проведения термодинамических расчетов. Результаты автора уже используются в некоторых исследовательских центрах (ИГЕМ РАН, ГИН РАН, ИЗК С-ПбГУ). а также широко внедрены в учебный процесс ЮТА.

Исследования, на основе которых подготовлена диссертация были сделаны автором либо самостоятельно, либо в составе возглавляемых им творческих коллективов. Исключение составляют работы по согласованию термодинамических свойств в системах As-S-0-Н и Sb-S-0-Н: здесь автор отвечал за теоретическую часть проблемы - разработку методов решения обратных термодинамических задач и построение на этой основе физико-химических моделей с участием компонентов системы.

Апробация работы. Результаты работы обсуждались на Всесоюзных конференциях Термодинамика в геологии (Миасс. 1988; Новосибирск, 1992). II Всесоюзном совещании Физико-химическое моделирование в геохимии и петрологии на ЭВМ (Иркутск, 1988). Всесоюзном семинаре Применение ЭВМ при гидрогеохимическом моделировании (Ленинград. 1991), международной конференции Thermodynamics of Natural Processes (Новосибирск, 1992). международной научной конференции Геофизика и современный мир (Москва. 1993), XIII международной конференции ИЮПАК по Химической термодинамике (Клермон-Ферран, Франция. 1994), XVI конференции Международной Минералогической Ассоциации (Пиза. Италия, 1994), и ежегодных конференциях МГГА Новое в науках о Земле (1986-1994). Содержание диссертации отражено в 28 публикация:-:.

Автор считает сьо;;м долгом выразить глубокую благодарность А. В. Зотову, оказавшему неоценимую помощь на стадиях постановки задач и

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

Диссертация состоит из введения, семи глав, заключения и дополнений общим объемом 363 страниц, включая 40 таблиц. 42 рисунка и список литературы из 196 названий. В дополнения включены' тексты расчетных программ.

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

Упрощенный вариант алгоритма положен в основу программы ■ SOL. предназначенной для компьютерной обработки результатов анализов природных и техногенных вод.

В настоящее время существует два подхода к расчету равновесий в многокомпонентных системах. Первый основан на том, что равновесию системы отвечает равновесие всех возможных химических реакций в системе AfG! = 0. Тем самым, этот метод основан на выполнении уравнений действующих масс или. иначе говоря, уравнений констант равновесия.

Альтернативный метод расчета основан на прямой минимизации функции Гиббса G системы для заданных давлении и температуре. Метод прямой минимизации не требует перебора возможных реакций и поэтому является значительно более "экономным" и гибким, чем метод констант, особенно для многокомпонентных систем. В диссертации для решения задачи расчета гетерогенных равновесий использован именно этот подход.

Минимизация функции Гиббса системы G, состоящей из J (Kj<J) компонентов водного раствора и Р (2<р<Р) "чистых" фаз

J X! Р

G = I xWjb + RTln—- + RTlriK,) + I Мп-Хр J-l X, p=2

при ограничениях типа баланса масс I (1<1<1) независимых компонентов J Р

I a13Xj + I at ,(j»p-j)-Хр = bj

j—i p-^

проводится с использованием метода неопределённых множителей Лагран-жа. В этих уравнениях р - стандартные химические потенциалы компонентов системы при заданных давлении и температуре Т [К]. К3 - коэффициент активности компонента j. х3 и Хр - количества [моль] компонента раствора j или фазы р в системе. Элемент матрицы ап отражает стехиометрию компонента 1 по независимому i. а b( [моль] - общее количество независимого компонента 1 в системе, равное его начальному значению.

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

В описание алгоритма включены приемы, позволяющие учесть окислительно-восстановительные реакции, ввести поправки на неидеальность компонентов раствора и проводить расчеты с вполне подвижными компонентами. Приведено обобщение алгоритма для расчета равновесий в мультисистемах. когда число возможных фаз системы превышает число независимых компонентов (системы с частично известным фазовым составом) . В последней версии к программе подсоединен широко известный пакет SUPCRT92. так что процедура расчета термодинамических свойств компонентов в широком интервале температур (до 1000°С) и давлений (до 5 кбар) проводится непосредственно в течение сеанса работы с BALANCE.

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

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

Во второй главе предложен и обоснован подход к описанию кинетики растворения минералов.

В настоящее время экспериментально наиболее хорошо изучены процессы взаимодействия силикатов и алюмосиликатов с водными растворами. а наблюдаемую кинетику растворения удается аппроксимировать зависимостью С(и ~ \}/г (параболический закон) [Ьаваяа А. С.. 1984]. Попытка объяснить параболическую стадию растворения явилась причиной появления диффузионной модели [Юо11аз1 И.. 1967]. Согласно этой гипотезе внешний слой минерального зерна частично изменяется так. что лимитирующей стадией становится диффузия вещества в объём раствора через этот изменённый ("выщелоченный") слой. Однако тщательное изучение поверхности минералов после их взаимодействия с растворами показало, что предполагаемые сторонниками диффузионной гипотезы вторичные слои, по-видимому, не существуют.

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

На рис. 1 изображены потенциальные барьеры растворения (44) и конденсации (-М). Для перехода иона в раствор необходимо затратить энергию, превышающую энергию связи иона в решётке. Константа скорости растворения определяется энергией активации Гиббса ДС~(«-М) и выражается уравнением Эйринга

ДС (-М)

К(—М> ~ ехр{--) , (2.1)

ИТ

г

где индекс """ относитсл к переходному состоянию активированного комплекса. Для осуществления обратного процесса ион должен сбросить

Рис. 1 Профиль функции Гиббса процессов растворения-"конденсации" на границе раздела фаз.

сольватирующие молекулы и промигрировать через двойной электрический слой

ДС (-М)

к(-М) ~ ехр{--} (2.2)

!?Т

Направление стрелок (-М) в формулах (2.1) и (2.2) показывает движение иона по отношению к межфазной границе М.

Рассмотрим минерал состава АХВУ. в котором энергии связей ионов в кристаллической решетке различны. В начале процесса растворения слабосвязанный ион покидает решетку быстрее и. таким образом, поверхность твердой фазы приобретает некоторый заряд, что эквивалентно появлению разности потенциалов на границе двойного слоя:

ДФ = Ф(М) - Ф(3) .

где Ф(М), Ф(Б) - потенциалы Гальвани поверхности и растворителя. Появление разности потенциалов ДФ является причиной изменения функций Гиббса и. следовательно, констант скоростей прямого и обратного процессов. Если Ф(М) < Ф(3), то вклад в АС (->М) отрицателен: катионы имеют тенденцию двигаться в область более отрицательного ■ потенциала. Принимая для простоты модель Гельмгольца двойного электрического слоя, получим, что функция Гиббса активации в случае движения к поверхности определится как

Айц"(-»М) = с" (-М) + а-р-а-ДФ . где Г - число Фарадея, 0 - заряд иона, а а ~ 0.5 - фактор симметрии, определяющий положение переходного состояния. Новая константа "конденсации" теперь примет вид

ДСа"(-М)

М-М) ~ ехр{--}

1?Т

или с использованием (2.2)

ц-Р-аДФ

кл(-М) = к(~»М) -ехр{ - - } . (2.3)

ИТ

а скорость

Б

Уа(->М) =--С-Ко(-М) (2.4)

V

Здесь Б - площадь поверхности минерала, V - обьём раствора, а С -концентрация частиц в растворе. Аналогично для процесса растворения "новая" константа скорости

О-Г- (1-а) -ДФ

ка(-М) = к(-М) -ехр{ - - ) . (2.5)

ИТ

и скорость

Б

У0(-М) ---С-кц(-М) (2.6)

V

Выражения (2.3-2.6) позволяют получить полные кинетические уравнения для процесса растворения минерала АХВУ. Вводя для удобства обобщённые безразмерные концентрации ионов А и В

. в

- Г

кА(-М) -кв (-М) кА (-М) -кв (*-М) -1

1/2

. в

Э Р

время т ---[к4 (-М) кп (-»М) ]'/г ■ С и потенциал <р = ДФ--и обозна-

V ИТ

чая и = |-1 , Г) = |-1 . получим окончательно

кв (-М)-1 квС-М)-1

(Зт

= - и-зеА-ехр{-алал1р) + Г)ехр{0А(1-аА)ф)

бэев 1 1

— ----эев ехр{-Овав«р} +--ехрШв(1-ав)!(>)

бт и Л

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

с1т

+ 0В

Дзев йт

= О

(2.7)

В диссертации путем численного решения системы (2.7) рассмотрены различные частные случаи растворения. Так. при 0А = йв, эеА = эев = зе в течение всего процесса растворения, и обобщенную концентрацию эе можно интерпретировать как степень протекания реакции растворения минерала (0<эе<1 и эе=1 при равновесии). Решение (2.7) в этом случае даёт

йзе <3т

1

эг

(2.8)

[(зе + 1/1)П) • (ае+ЦП)]1/г

Параметром, определяющим кинетику накопления частиц, является фактор ЦП, который характеризует "асимметрию" поведения ионов А и В в кристалле и растворе. При ЦП = 1 йэеЛЗт = 1 - эе и кинетика описывается обычным экспоненциальным законом. Однако при нарушении "симметрии" начинают проявляться некоторые особенности: а) в начале растворения (зе < 0.2) наблюдается резкое уменьшение скорости, что в рамках нашей модели связано с возрастанием потенциала на границе раздела фаз: 61 кинетика накопления частиц на начальном участке (эе < 0.5) достаточно хорошо аппроксимируется параболическим законом; в) при приближении к

равновесию закон скорости приобретает линейный характер.

В рамках предложенной модели получили количественное объяснение экспериментальные факты зависимости скорости растворения силикатов (Мя25104, форстерит и М^103> энстатит) от кислотности раствора, а также предложено решение широко^ известного вопроса об образовании арагонита в морских осадках.

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

В третьей главе предложен метод оценки электрострикционной составляющей химического потенциала иона в водном растворе. Метод основан на сравнении энергий взаимодействия диполей Н20 до и после введения иона в раствор. Полученное на основе этого подхода уравнение состояния иона в растворе совместно с электростатической составляющей модели Ш" позволяет добиться удовлетворительного согласия с экспериментом в широком интервале давлений и температур.

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

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

- М°0 + М°1»1 •

где - "собственный" химический потенциал компонента, а -

член, характеризующий энергию взаимодействия в системе ион-раствор. Свободная энергия этого взаимодействия включает в себя работу заряжения иона в растворе (/1°е - электростатическая часть) и изменение энергии растворителя при Бнесении в него иона (м°„)

М°1П1 = М°п + М°е

- И -

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

й2 • е2 • N. 11 й2 1

мое -----(_ _ 1)ш л---(_ _ 1}

8Я-£0 Ге £ Ге I

где а - заряд иона, с - диэлектрическая проницаемость чистого растворителя, ге [А] - эффективный электростатический радиус иона в растворителе, а п = е2-На/8яе0 = 1.66027-Ю5 кал-А/моль. Вклад р.°е в становится определяющим при повышенных температурах и невысоких давлениях, когда величина с невелика. Вводя вслед за авторами модели НОТ обозначение для коэффициента Борна частицы

О2•е2• 1 й2 со -- • - = П" -г- .

8я-Е„ ге геШ

получим

1

м°е - Ш-( - -1) (3.1)

£

Величину парциального мольного объёма компонента можно получить дифференцированием его химического потенциала по давлению Чг = (с1ц0/<1Р)т, так что

Vг - V. + V. + Уе .

где

Уе = (йц°е/йР)т = -и-О + (1/Е - 1)(бш/бР)т электростатическая составляющая объёма. 0= -3/£2 ■ (йе/(ЗР)т, у0 = (а/1°0/йР)т - "собственный" объём частицы. V,, = ((Зд°п/бР)т - изменения в объеме растворителя, вносимые ионом (электрострикционная составляющая).

При низких температурах вклад в величину Уг является заметным и именно он определяет температурную и барическую зависимости объёма частицы.

Для оценки V,, привлекались различные подходы. Все они однако содержат физически мало обоснованные подгоночные параметры. Так. ещё со времён Таммана контракцию растворителя при внесении в него иона связывали с наложением некоторого "внутреннего" давления Ре. В этом случае (в расчёте на 1 кг растворителя)

У„(Р) = ГVСР+Ре) - У(Р)3 -55.51 Используя уравнение Гэйта для жидкой воды V(Р) = А - С-1п(В+Р).

где У(Р) - мольный объем НгО при давлении Р; А. В и С - постоянные не зависящие от давления, легко получить для не очень высоких Ре

У„(Р) « - Г - 1п(---) = С'-—

В+Р В+Р

Здесь С" = С-55.51, С = 2.74 см3/моль. В = 3525 бар для Т=25°С.

Это уравнение лежит в основе практически всех количественных методов предсказания барической зависимости мольных объемов электролитов Шикег^ее Р.. 1961; Львов С.Н.. 1983; Любимов С.Л.. 1985]. Выбором Ре для каждого электролита удавалось добиться хорошего описания экспериментальных зависимостей. Тем не менее, величина Ре является чисто подгоночным параметром. Кроме того, использование такого феноменологического подхода не позволяет построить правдоподобной модели для температурной зависимости внутреннего давления Ре.

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

В диссертации рассчитаны энергии диполь-дипольного взаимодействия в растворителе до и после внесения в него иона. Общее изменение энергии взаимодействия 1 моль диполей Нг0 при введении ионов в раствор составит

Ди„ = Е2 • - . —1 - н - — (3. 2)

зщ г6 г3

Здесь Н = - = 2.038-Ю5 Дж-А3/моль. п2 и п, - число молекул Нг0,

4Я£0

окружающих данную до и после внесения иона, а г имеет смысл среднего расстояния между взаимодействующими диполями. Величину г естественно связать с мольным объёмом "невозмущённого" растворителя V ~4/Зл-иА - (г/2)2, откуда (г [А])3 -3.172 -V [см3-моль"1 ]. Подставляя в (3.2), получим

п, п,

Дит [Дж/моль] - 6.42-104 --- 1.65-108 • - (3.3)

V ТУ2

Общее изменение свободной энергии при внесении 1 моль ионов в растворитель тогда

М°„ - гп-Лит . (3.4)

где т - число молекул Н20. находящихся под воздействием одного внесенного иона. Можно показать, что "действие" иона ограничивается первой гидратной сферой, а ш. таким образом, должно коррелировать с гидратным числом иона.

Теперь, подставляя (3.3) в (3.4). получим

6.42-101 -п2 1.65-108-П1

- ш ( -- - --- ) -

V Т-У2

1 2. 57-103 -П, /пг = т-6.42-104---(1 - -) .

V Т-У

гдь второе слагаемое в скобках характеризует изменения "чисел окружения" взаимодействующего диполя при внесении в растворитель иона.

Изменение объема растворителя в расчете на 1 моль ионов теперь можно найти дифференцированием по давлению

6. 42-10* -По ¿IV 2-2.57-103-П,/Пг УЛсм3/моль]~ -10-го---(—)| [1--]

V2 6Р |т ТУ

Для расчёта объема "невозмущенного" растворителя предлагается использовать уравнение Тэйта в виде, пригодном для широкого диапазона температур V = А - С-1п(В0 + Р + Вг-Р2) [Ходаковский И.Л.. 1975]. поэтому

(IV 1 + 2 • В2 ■ Р

(—) I = С

ар 1Т в0 + р + в2-р2

6.42-пг-С 1+2 Вг-Р 5. 14-103-П./пг

Уп ~ 10-Щ---------[1---——] (3.5)

V2 В0 + Р + В2-Рг У-Т

Анализ уравнения (3.5) показал, что в диапазоне экспериментальных условий 0 - 100°С и 1 - 5000 бар. где вклад \'п значителен, зависимости С(Т) и V(Т. Р) слабо влияют на величину \'п. так что полагая их постоянный; запишем окончательно

1+2 • В, Р а, Ю т а!---[1--] . (3.6)

В„ + Р + В2-Рг Т

где а, и а2 - некоторые характеристики растворителя, не зависящие от температуры и давления, а В0 и В2 - параметры модифицированного уравнения Тэйта. зависящие только от температуры.

Уравнение (3.6) сохраняет характер аддитивности при . переходе от иона в растворе к электролиту, так что в дальнейшем оно было использовано для "настройки" модели: нахождения параметров а, и а2 с использованием экспериментальных значений для ИаС1 ( а! = 0.9669-102 кал/моль и а2 = 4.902-10г К). Результаты расчёта для серии электролитов представлены на рис. 2 и позволяют сделать вывод о хорошем согласии с экспериментом. При этом "электростатические части" объемов рассчитывались с использованием значений ш по модели НЮ", а величины У0 и ш выбирались таким образом, чтобы удовлетворить эксперименту при Р = 1 бар для двух температурных точек: 0 и 25°С.

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

1 + 2 - в2 - Р а2 1 с1ш У2 = У0 + Ю-т-а,---[1--] - ш-а + (- -1) • (—) I .

В0+ Р + В2-Р Т е С1Р |т

пригодное для широкого диапазона давлений и температур. Уравнение содержит всего два подгоночных параметра. V,, и т. каждый из которых имеет ясный физический смысл: У0 - "собственный" объём иона, который может быть независимо оценен из данных по плотности расплава соли; ш - число молекул Н20. входящих в область "влияния" иона. Значение га может быть найдено из единственного экспериментального значения У2 иона при стандартных условиях. Использование такого подхода может оказаться весьма полезным для предсказания термодинамических характеристик ионов при повышенных температурах и давлениях в условиях дефицита экспериментальных данных об этой частице.

В четвёртой главе предложено описание процессов ассоциации 1:1 электролитов, не требующее индивидуальных подгоночных параметров.

У 2 смЗ/моль

кш,

V 2 СМЗ/МОЛЬ

Ка2804

а) Т°С

Р бар

Рис. 2. Температурная (а) и барическая (б) зависимости парциальных мольных объемов некоторых электролитов. п-.'.к;'» - расчёт по уравнениям (3.2), (3.7), точки - эксперимент

Подход включает в себя алгоритмы расчета неэлектростатических и

сольватных свойств ионной пары. Получены выражения для параметра Борна ассоциата. состоящего из двух неполяризованных гидратированных ионов. Модель "состыкована" с уравнением состояния НКГ и поэтому пригодна для расчета любых термодинамических свойств ионной пары.

В настоящее время для интерпретации экспериментальных данных по константам устойчивости ионных пар используют практически два различных подхода. В рамках "плотностной" модели [МагзПаП IV. Ь.. 1970. 1972. 1991] считается, что константа равновесия реакции ассоциации определяется плотностью растворителя. Однако физические обоснования такого предположения нельзя считать достаточно убедительными, и поэтому даже сегодня эту модель следует рассматривать как полностью эмпирическую.

Другой подход носит название электростатического и восходит к работам Бьеррума [ВЗеггит N., 1926]. Согласно этой модели взаимодействие противоионов определяется исключительно кулоновскими силами притяжения. При этом считается, что оба иона, образующие ионную пару все еще окружены молекулами растворителя и каждый из них сохраняет свою электростатическую энергию. В этом случае функция Гиббса реакции ассоциации двух ионов А и В с зарядами йА и С1в эквивалентна энергии кулоновского взаимодействия ионов, находящихся на расстоянии а в среде с диэлектрической проницаемостью е:

где Na = 6.022045-Ю23 моль"1 - число Авогадро. а ц ---

8Xt 0

= 1.60027 -Ю5 кал -А/моль. Это простое выражение и представляет собой основное уравнение электростатической теории [FuossR.M.. 1958]. Дальнейшее его усовершенствование возможно путем учёта сил отталкивания в ассоциате [Мелвин-Хьюз Е.А.. 1972], учёта влияния плотности растворителя на энергию взаимодействия ионов [Рыженко Б.Н.. 1981. с. 118]. учета поляризации ионов [Брызгалин 0. В.. Рыженко Б.Н.. 1982] или путем ввода дополнительных членов, характеризующих неэлектростатическое взаимодействие [Brady P.V.. 1990]. Важно отметить однако, что общим для всех подходов является эмпирический характер межионно-

AG =

NA-QA-QB-e* J_ 4ji-E 0е а

(4.1)

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

Вслед за подходом HKF любую термодинамическую характеристику водной частицы представим в виде суммы

= = Н„ + Не .

где Ее - "сольватная" составляющая, которая отражает вклад электростатического взаимодействия между ионом и растворителем; Sn - "структурная" составляющая, которая отражает вклад остальных, в основном ковалентных взаимодействий в системе частица-растворитель.

Общей чертой практически всех электростатических моделей является допущение о независимости неэлектростатической части функции Гиббса комплекса от температуры и давления. Это означает, что неэлектростатические составляющие энтропии и теплоемкости ассоциации принимаются равными нулю для всех Т и Р: ArSn = 0; ArCP п= 0. Предлагаемый в данной работе подход позволяет отказаться от этого допущения и фактически представляет собой метод оценки ArS„ и ДгСРп. Тем самым развиваемый подход можно рассматривать как обобщение электростатической теории.

Для расчета структурных составляющеих реакции ассоциации ионов А + В = АВ (4.2)

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

G° = G°(0) - RT(I ln q(m) - ln NA). где T - температура [К], G°(0) - функция Гиббса частицы при Т=0. q(m) - функция распределения молекулярной моды ш , а суммирование ведется по всем модам молекулы.

Поскольку S = -(dG/dT)|p . получим

1 da,m)

S = Ш ln q(ra) -ln N.) + RT-lnd ----:-)

q(m) cT

Последний член выражения может быть преобразован с учетом того, что

он имеет смысл среднего запаса энергии в молекуле. Поэтому

1

S = SCT. + S3H. = R(I ln q(m) -ln Na) + —-ZNa<E>.

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

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

для поступательного движения частицы массы М [г/моль]:

q(t) = 1.88-1023-Т3/г-М3/г. (4.4)

для вращения линейной молекулы АВ вокруг центра масс:

qAB(r) = 4.123-10~2- ( * ) t2-T . (4.5)

где I [А] - длина связи А - В; АВ для осциллятора с частотой v [Гц]:

qAB(v) - --(4-6)

1 - ехр{-а}

где а = 4.798-Ю"11-v/T - характеризует "жёсткость" связи, в диссертации получено соотношение для структурной составляющей энтропии реакции ассоциации

МАВ3/г 1

AS„ = R- 1п[7.65-10~3 ■ 1г------] +

МА3/г-Мв3/г 1-ехр{-а)

ехр(-а)

+ R- [а----0.5] (4.7)

l-exp(-a)

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

ехр{-а/2)

ДС„ = R- [а---0.5]

1-ехр(-а)

(4.8)

Считая, как и прежде, что образование, ионной пары представляет собой простое объединение двух гидратированных фрагментов, будем считать, что ДУП = О. (4.9)

Как уже отмечалось выше, при описании сольватной части химического потенциала авторы модели НКГ используют уравнение Борна (3.1). которое имеет смысл энергии гидратации иона с зарядом 0 и радиусом ге в среде с диэлектрической проницаемостью е, а ш = ч1-а2/ге - параметр Борна этого иона. Важно отметить, что под £ здесь понимается диэлектрическая проницаемость чистого растворителя (Н20) в глубине раствора. Таким образом, все поправки, связанные с "континуальностью" модели и эффектами диэлектрического насыщения вблизи иона сведены в величине ге. которая интерпретируется как электростатический радиус иона, и по существу должна пониматься как эмпирический параметр. В рамках модели НКГ однако установлены жесткие зависимости мевду кристаллографическим и электростатическим радиусами, так что для "свободных" частиц эта величина оказывается определенной и не является индивидуальным подгоночным параметром.

В последнем пересмотренном уравнении НКГ (1988) для улучшения точности описания учитывается зависимость электростатического радиуса от давления и температуры согласно уравнению

ге - гх + 101 -{к + g(P.T)} Зависимость от Р и Т определяется корректировочным множителем Я(Р. Т). который одинаков для всех водных частиц и также не является индивидуальным подгоночным параметром.

Выражение (3.1) как бы эквивалентно включению электростатической теории в уравнение состояния водных частиц. Например, можно показать. что соотношению (4.1) соответствует изменение параметра Борна в результате образования ассоциата с эффективным межионным расстоянием а

Оа-ОВ

Ды = 2ц----(4.10)

а

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

ностью не худшей, чем экспериментальная ошибка, однако его предсказательная способность, видимо, ограничена. Так. принятые значения параметров Борна многих ассоциатов 1:1 электролитов к = 0. что соответствует физически мало реальной ситуации близкого к нулю дипольно-го момента ионной пары.

В диссертации предложен метод расчета параметра Борна ассоциата (оАВ. состоящего из двух неполяризованных сферических ионов. При этом автоматически учитывается уменьшение напряженности электрического поля вблизи иона, вносимое за счет присутствия противоиона. Полученное точное выражение шАВ ассоциата достаточно громоздко и поэтому здесь не приводится. Для нас интерес представляет величина Дш = шАВ - шА - шв изменения параметра Борна в реакции ассоциации (4.2). В работе показано, что значение Дш можно аппроксимировать уравнением (4.10). где, однако, а ~ 0.7338-1 при малых значениях межионного расстояния, а при больших I а ~ 0.5-1. Этот результат объясняет низкие значения параметра а. получаемые обычно при подгонке к эксперименту.

Уравнения (4.7-4.9) использовались нами в дальнейшем для описания неэлектростатических термодинамических свойств реакций ассоциации. Для расчета электростатической составляющей функции Гиббса реакции использовали выражение аналогичное (3.1), в котором Дше рассчитывали по точному варианту уравнения (4.10). При низких температурах межьядерное расстояние принимали равным сумме электростатических радиусов компонентов ионной пары I = ге(А) + ге(В), а при повышенных параметрах состояния в эти значения вводили поправки. Проведенный анализ показал, что возрастание (по абсолютному значению) величин Дшр т в малоплотных флюидах может определяться уменьшением межатомного расстояния I в ионной паре. .Действительно, как видно из уравнения (4.Ю) "стягивание" ионной пары, обусловленное перекрыванием гидратных оболочек, приводит к возрастанию |Дш|. причем в диапазоне межионных расстояний 1е = ге(А)+ге(В) и 1Х = гх(А)+гх(В) эта зависимость почти линейна. В таком случае для параметра Борна реакции ассоциации ионной пары До) с межионным расстоянием I можем записать

I -1е

Аш( I) = Дше +■ (Д- Дше)--. (4.11)

1х-1е

где Дсое и Дшх - параметры Борна. соответствующие случаям "соприкаса-

ющихся" (1е) и "плотно сжатых" (1х) заряженных сфер.

Предположение о "стягивании" ионной пары естественно связать с дегидратацией ионов в малоплотных флюидах: при уменьшении плотности растворителя равновесие реакции

А-(Н20)п = А-(Нг0)„.] + Н20 будет сдвигаться в сторону продуктов. Уменьшение гидратного числа иона в этом случае будет способствовать образованию более "тесных" ассоциатов, что в свою очередь эквивалентно уменьшению среднего межатомного расстояния в ионной паре.

Для описания зависимости Д1 = I - 1е от плотности растворителя в диссертации предложено использовать корректировочный я-фактор модели ШТ. В этом случае предлагаемая нами модель достаточно естественно встраивается в формализм Н№ и становится пригодна для расчета любых термодинамических свойств ассоциатов. Проведенный анализ показал. что для Д1 можно использовать соотношение

Д1 = 2g•Цx - 1е) , (4.12)

так что для расчета параметра Борна реакции согласно (4.13) будем иметь

Да) = Дй)е - 2g• (Дшх - Дше) (4. 13)

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

Таким образом, предложенная модель для описания ионной ассоциации сильных электролитов

- не содержит индивидуальных подгоночных параметров;

- "состыкована" с уравнением состояния НКГ и поэтому пригодна для расчета любых термодинамических свойств ионной пары;

- адекватно описывает поведение 1:1 электролитов в широком диапазоне сверхкритических параметров.

В пятой главе методом физико-химического моделирования на ЭВМ проанализированы термодинамические свойства компонентов "типичных" гидротермальных систем Аэ'11-Б11-0-Н и БЫ 1'-Б11-С1-0-Н в интервале температур 25-350°С с использованием уравнения состояния НКГ для частиц в водном растворе. Получены согласованные значения термодинамических параметров веществ, отличающиеся от имеющихся в сводках.

Таблица 1. Рассчитанные значения термодинамических свойств реакций ассоциации различных 1:1 электролитов: ДСРГ Тг [кал/моль] - стандартная функция Гиббса при Рг = 1 бар и Тг = 298.15 К: ДБ,, [кал/(моль-К)], ДС„ [кал/(моль К)] - неэлектростатические составляющие энтропии и теплоемкости при стандартных условиях; До)е-10"5, Дшх-10"5 [кал/моль] - параметры Борна для межъядерных расстояний равных суммам электростатических ге(А)+ге(В) и кристаллографических гх(А)+гх(В) радиусов ионов, составляющих ассоциат.

Электролит Аврг.тг АБ„ ДС„ Дше•10"5 Дшх • 10" Б

ЫаС1 560 -4.704 0.973 -1.2162 -1.3917

ИаВг 1080 -4.435 0.979 -1.1219 -1.2832

Ш 830 -4.020 0.982 -0.9921 -1.1440

КС1 1265 -4.071 0.982 -1.1098 -1.2419

КВг 1128 -3.824 0.988 -1.0697 -1.1903

К1 1458 -3.442 0.988 -1.0118 -1.1183

СэС1 -452 -3.520 0.987 -1.0265 -1.1448

СэВг -158 -3.290 0.990 -0.9915 -1.0985

Сэ1 -122 -2.936 0.992 -0.9408 -1.0342

-1270 -4.669 0.975 -1.2211 -1.5200

ИэС1 47 -3.840 0.986 -1.0737 -1.1995

ИЬВг 504 -3.601 0.989 -1.0359 -1.1502

ИМ 730 -3.230 0.991 -0.9812 -1.0817

ЫС1 1717 -5.226 0.931 -1.2631 -1.5026

Результаты могут быть использованы для описания поведения систем в широком диапазоне температур (до 500°С) и давлений (до 3 кбар).

В литературе имелись довольно разрозненные сведения по поведению Аб1П и БЬ111 в водных растворах, а разные термодинамические сводки давали весьма различные значения свойств Аэ- и БЬ-содержащих компонентов. Кроме того, даже стехиометрия некоторых комплексов не была установлена. В такой ситуации наша задача состояла в выборе наиболее надежных и внутренне согласованных экспериментальных данных, их новой интерпретации (если это было необходимо), а в некоторых случаях - в проведении дополнительных уточняющих экспериментов по растворимости (все эксперименты проводились Зотовым А. В. и сотрудниками лаборатории физико-химического моделирования ИГЕМ РАН). Полученные результаты были использованы для восстановления термодинамических параметров конкретных частиц и согласовании этих значений.

3 5

О

6 ■

ТС

400 600 800 СбВг

ГС

400 600 800

400 600 800

Рис. 3. Экспериментальные (кружки) и рассчитанные по уравнениям (4.8-4.13) (линии) значения 1вКаз5 реакций ассоциации различных электролитов в зависимости от температуры при разных давлениях. Размер кружка соответствует экспериментальной ошибке. Цифры на рисунках отвечают давлениям в кбар. а Рнас

- давлению насыщенного пара н20. табл. 1.

Данные для расчета взяты из

Процесс согласования представлял собой итерационную процедуру, на каждом шаге которой рассчитанные значения растворимостей сравнивались с "эталонными" экспериментальными. Первоначальные оценки для параметров уравнения состояния водных частиц получали с использованием развитых в работе методов. Эти значения служили основой для расчета термодинамических свойств при повышенных температурах и давлениях. Расчет растворимостей проводился по программе BALANCE, a HKF параметры водных частиц восстанавливались с использованием программ точной оптимизации UT-HEL (автор Шваров Ю.В.) и авторской HELGES0N.

Имевшийся экспериментальный материал был, как правило, ограничен невысокими значениями температур (до 200°С) и низкими давлениями насыщенного пара Н20. Этих данных явно недостаточно для моделирования гидротермальных процессов с участием As- и Sb-содержащих частиц. Таким образом, остро стояла задача корректной экстраполяции их термодинамических свойств на область более высоких параметров состояния (до 500°С и 3 кбар). Именно поэтому в диссертации термодинамические параметры водных частиц определялись с использованием уравнения HKF. Такого рода сводка термодинамических свойств частиц является первой для систем As111- SI!-0-H и Sb^-S^-O-H (табл. 2, 3).

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

Типичной задачей геохимического исследования является анализ форм нахождения и миграции минерального компонента в гидротермальном растворе, а также выяснение физико-химических условий формирования оруденения или преобразования рудных минералов. В шестой главе на основании полученных термодинамических данных обсуждается геохимическое поведение As и Sb в гидротермальных растворах. Построены диаграммы полей преобладания различных As- и Sb-содержащих частиц в зависимости от рН. концентрации сульфидной серы (mZH2S) и летучести водорода (fH2) при 25 - 300°С.

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

Табл. 2. Стандартные термодинамические свойства соединений в системах Азп1.3п_0_н и 5Ьп1.5п_о_н

Вещество Г ,298 кал/моль со ь 298 кал/моль/К У°2 98 см3/моль Ср(кал/моль/К) = а+Ь-Т-с/Тг+й-Т2 Температ. диапазон °С

а Ь-103 с•10"5 сЫО6

Аэ(кр.) 0 8.5206 12.963 5.6341 0.7125 0.400 0 25-727

Азг03(мон.) -138009 29.2854 42.6530 -27.0507 8.6501 18.309 25-314

АЭгОз(Куб.) -137837 25.6644 51.32 6.6635 75.0478 0 -66.1807 25-278

АзгБ3(аури) -21866 39.0535 70.52 42.7608 -52.6250 -4.3098 60.2861 25-312

Аз4Б4(реа) -31721 58.68 119.2 51.77 -0.053 5.58 0 25-903

5Ь(ромб. ) 0 10. 773 18.178 6. 6876 -2. 9689 0. 1659 4.8506 25-267

5Ьг03(сен) -149967 30.091 52.21 28.2605 -1.9551 3.0471 8.4419 25-879

5Ьг0з(Вал) -148193 32.118 50.01 30.00 -1.9551 3.0471 8.4419 25-928

БЬ^з (ант) -34156 42.854 73.414 23.889 13.2526 -0.7192 0 25-822

О!

Для ' и " полное уравнение теплоемкости имеет вид Ср = а + Ь-Т - с/Т2 + й'Т2 + е-Т3, где

' (МО6 = 2. 1713-Ю"3, е-Ю10 - -3.3137 И " б'Ю6 - 2.7759, е-Ю10 - 6.8024

Табл. 3. Стандартные термодинамические свойства водных частиц в системах Asn i _3П -о-Н и Sb111-S11 -0-Н

Соединение ¿G°r,298 S°298 ш-iO5 V°298 r 0 298 a!-10 а2-Ю-2 аз a4-10"4 Cj <V10-"

a 6 a B 6 Г a Д e 6 ' e

H2As03" -140330 26.4 1.2305 26.02 -2.11 5.7934 6.3646 3.2485 -3.0421 15.8032 -3.6253

H3As03° -153005 47.80 -0.4111 34.52 29.06 6.3337 7.6759 2.7294 -3.0966 19.4362 2.8931

HAs038" -123788 -3. 600 3.2971 8.78 -63.13 4.2208 2.5215 4. 7586 -2. 8824 -3. 7835 -14.422

As033" -105501 -44.598 5.5349 -19.5 -113.0 1.1991 -4.8540 7.6580 -2.5779 -9.3841 -26.154

AsaS3° -13066 24.300 -0.0301 28.51 34.60 5.6549 6.0261 3.3810 -3.0282 26.1711 4.0134

HAS2S4" -16995 52. 008 0.8413 51.68 31.04 9.1539 14.6033 0.0115 -3.3939 32.0507 3.2744

ASzS42- -5711 53.012 2.4140 65.68 -27.89 11.6635 20.7266 -2.3901 -3.6329 11.87.86 -8.7476

Asi51¿" -24530 104.995 1.6228 116.0 12.32 18.2600 36.7830 -8.7237 -4.3021 28.3461 -0.6095

Sb(0H)3° -153890 47.0 -0.23 42. 18 39.49 7.4932 10.3149 1.6163 -3.2137 27.207 5.0115

Sb2(0H)6° -307885 104.1 -0.46 84.73 79.00 14.9864 21.0298 3.2326 -6.5462 54.4414 10.023

Sb(OH)4" -194320 48.0 0. 8946 48.00 0.96 8.6736 13.397 0.4833 -3. 3328 14.8793 -2.839

Sb2S42" -14000 46.0 1.8153 47.34 -51.80 8.9334 14.0315 0.2339 -3.3591 -7.6536 -13.586

HSb2S4- -27350 63.33 0.7274 63.74 -82.70 10.7077 18.3638 -1.4692 -3.3582 -35.674 -19.881

H2Sb2S4° -33800 82.70 -0.9124 74.20 -11.60 11.5711 20.4729 -2. 2984 -3.6254 -8.944 -5.3975

Примечание: а-кал/моль; б-кал/моль/К; в-см3/моль; г-кал/моль/бар; д-кал-К/моль/бар; е-кал-К/моль

I 4

ре. Однако они скажутся на условиях кристаллизации сульфидных фаз, так как понижение окислительно-восстановительного ютенциала приводит к тому, что вместо аурипигмента устойчивыми фазами становятся реальгар (Аэ^) и самородный мышьяк. Построенные диаграммы позволяют проследить относительную устойчивость этих фаз б зависимости от Гнг. рН. т1Н23 и температуры. Диапазон изменения лпучести водорода при этом был выбран с учётом реальных условий на с:временном мышь-яково-сурьмяном рудопроявлении в кальдере Узон на Камчатке.

Для ЭЬ111-Б11-С1-0-Н системы показано, что с р:стом температуры всё большее значение приобретают гидроксокомплексы, а поле стабильности антимонита резко сокращается и смещается в более кислую область. Таким образом, при температуре выше 250°С основной формой нахождения сурьмы в гидротермальных растворах являега ЗЬ(0Н)3°. Тем самым подтверждается общепризнанная точка зрения [Сорокин В.И.. 1988] о превалирующей роли этой формы сурьмы в гидротермальных процессах.

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

Температура по разному влияет на растворимссть антимонита в сульфидных растворах (с образованием сульфидных комплексов) и в растворах с преобладанием гидроксокомплексов. В пегзом случае с понижением температуры растворимость изменяется сравнительно слабо, а во втором, напротив, уменьшается очень резко. В последнем- случае падение температуры (в сравнении с другими факторами: подкислением раствора, привносом в систему сероводорода) являетсл. по-видимому, важнейшим фактором рудоотложения.

В работе построена модель рудоотложения в калькре Узон на Камчатке. позволяющая объяснить наблюдаемую смену мине:альной зональности (антимонит - реальгар - аурипигмент) и согла:овать расчетные величины концентраций водных частиц с измеренными гт з!1и значениями. Модель представляет собой дополненную сульфаг:едукцией модель охлаждения с окислением [Бычков А.Ю.. 1991]. Главней Факторами ми-нералообразования при этом являются охлаждение раствора и изменение

его окислительно-восстановительного потенциала за счёт процессов сульфатредукции.

Последняя седьмая глава диссертации посвящена разработке нового физико-химического подхода к описанию эволюции флюидов, богатых летучими компонентами. Начало вскипания такого флюида обычно рассматривают как мощный спусковой механизм начала рудоотложения. Косвенным подтверждением такой точки зрения являются многочисленные экспериментальные данные по составу газово-жидких включений, генетически связанных с оруденением: отмечается, что интенсивная минерализация имеет место при эволюции первоначально богатых С02 гидротермальных растворов. Главными физико-химическими факторами рудоотложения обычно считают изменение кислотности или уменьшение концентрации лиган-дов при удалении летучих компонентов. При этом, как правило, модели не рассматривают изменение диэлектрической проницаемости среды при изменении состава флюида. Буквально в последние годы появилось несколько работ [Walter J.V., 1990; Shott J.. 1991, 1992; Колонии Г.P., 1993], в которых отмечается, что это обстоятельство может оказаться источником грубых ошибок или делаются оценки такого влияния, но на равновесие отдельных реакций.

В диссертации предложен подход, позволяющий последовательно учесть влияние растворённых неполярных газов на реакции комплексооб-разования во флюиде. Использование такого подхода позволяет не ограничиваться рассмотрением отдельных реакций, а исследовать равновесия в гидротермальных системах в целом. Проведённый в рамках этого метода расчет растворимости галенита PbS в модельном H20-C02-NaCl-H2S флюиде позволил выявить весьма необычные химические свойства таких растворов. Оказалось также, что изменение диэлектрической проницаемости смеси при вскипании газонасыщенного флюида следует рассматривать как мощный фактор рудоотложения.

Как уже отмечалось выше (гл. 3. 4) для реакции комплексообразо-вания типа

А + В = АВ , (6.1)

изменение Функции Гиббса можно записать в виде

drG° = ArG°n + ArG°e = ArG\ + ДгшАВ( 1 /£ - 1).

I I

где ДГС°П = /1°„ (АВ) - ц"п(А) - 11°„(В) - "структурная" составляющая реакции ассоциации, а Дгшдв = ы(АВ) - и(А) - ш(В) - изменение параметра Борна в результате реакции.

Если величина ш(АВ) неизвестна, значение ДгшАВ можно оценить с использованием электростатической теории или по модели ионной ассоциации (4.10):

„ ЧаОВ

Дго)Ав = 2ц -

а

Легко видеть, что поскольку для противоионов йА-йв<0, ДгыАВ < 0. Увеличение температуры от ?! до Т2 обычно сопровождается уменьшением диэлектрической проницаемости растворителя (ет1<еТ2). поэтому

1 1 1

Д(ДгС°)=ДгСе(Тг)-ДгСе(Т1) = Дгш- (---) = Дг10-Д(-) < О,

ет2 ет1 е

и равновесие реакции комплексообразования смещается в сторону продуктов. Именно этот механизм лежит в основе объяснения устойчивости ассоциатов в околокритической области даже для "традиционно" сильных электролитов.

В нашем случае добавляется ещё одна причина изменения диэлектрической проницаемости флюида. Растворение неполярных газов (С02. СН4 И2) в воде также будет приводить к снижению диэлектрической проницаемости смеси е и. следовательно, дополнительно сдвигать равновесие реакций ассоциации в сторону продуктов. Вклад этого эффекта в изменение состава рудообразующего раствора нам и предстоит оценить.

Для расчета диэлектрической проницаемости смеси можно воспользоваться уравнением Ландау, в котором аддитивным оказывается кубический корень из е

Е1'3 = 1Ф1Е11/3 =Ф„го-ЕН2о1/3+ Фсог ' £сог1/3+ ФСН4 ' ^сн*1/3 +.-Л6.2) 1

Здесь Ф, - объёмная доля компонента 1 во флюиде с диэлектрической проницаемостью £,. Можно показать, что поскольку е неполярных газов близки к единице, даже небольшие их концентрации могут приводить к заметному снижению диэлектрической проницаемости смеси е. Это изме-

нение t может оказывать заметное влияние на распределение рудообра-зующих элементов в растворе и. следовательно, на растворимость рудных минералов.

Для расчета поправки к химическому потенциалу компонента можно воспользоваться соотношением Борна •

1 1

Ajli = ы-(--:-) (6.3)

£ ен20

Здесь £цго " диэлектрическая проницаемость чистого растворителя при давлении Р и температуре Т, а £ - диэлектрическая проницаемость смеси при заданной концентрации СОг и тех же Т и Р.

Необходимо также учесть, что изменение диэлектрической проницаемости среды приведет к изменению коэффициентов активностей ионов в растворе.

Для проведения количественных расчетов была выбрана система PbII-Na+-Cl"-C02-SII-H20, компонентам которой отвечает надёжная термодинамическая информация. Задача заключалась в расчёте растворимости галенита PbS в условиях, моделирующих эволюцию водно-углекислот-ного H20-C02-NaCl флюида.

Параметры среды минералообразования и состав рудообразующего раствора обычно восстанавливают исходя из состава флюидных включений. Эти параметры сильно варьируют, однако можно считать, что наиболее "продуктивным" является диапазон температур 400 - 150°С. Давления. во флюиде также меняются в широком диапазоне от 300 бар до 3 кбар. Нами была взята промежуточная величина Р = 1 кбар. Концентрация NaCl во флюиде была принята ~ 1 моль/кг. что также можно считать типичной для реальных рудообразующих гидротерм.

Первоначальная концентрация С02 во флюиде была принята равной 9.49 моль/кг Нг0. Эта величина соответствует температуре гетерогени-зации Тг = 350°С. Таким образом, эволюция модельного флюида представляла собой на первой стадии просто остывание гомогенного раствора с постоянной концентрацией С02 в диапазоне температур 400 - 350°С. При 350°С начиналось "вскипание" и концентрация С02 в растворе уменьшалась при снижении температуры, двигаясь вдоль кривой предельной растворимости. Концентрация H2S во флюиде была принята равной 0.01 моль/кг.

Для того чтобы выделить влияние эффекта снижения е смеси, раст-

воримость PbS и полный равновесный состав раствора рассчитывались в двух вариантах. В первом "традиционном" варианте все термодинамические данные были подготовлены с использованием пакета SUPCRT92, при этом неявно диэлектрическая проницаемость флюида принималась равной еН20 при заданных Т и Р. Во втором - "точном" - для расчёта химических потенциалов всех водных компонентов согласно (6.2-6.4) были внесены необходимые поправки. Расчёт растворимости PbS проводился с помощью программы BALANCE.

Результаты расчёта позволяет утверждать: учёт диэлектрической проницаемости кардинально меняет всю "химию" компонентов раствора! Во-первых, из почти нейтрального (рН=5.39) раствор становится щелочным (рН=6.72). Во-вторых, сильно изменяется концентрация "свободного" иона С1~ в растворе (здесь и дальше все концентрации и растворимости даны в моль/кг Н20): [СГ] (традиционный) = 0.7182; (СП (точный ) = 1.765 1О"2. Причиной этого является сильная тенденция к ассоциации для реакций

Н* + ОН" = Н20 Na+ + СГ ~ NaCl0 НСО3- + Н* = С02° + Н20 Поведение свинца в растворе также сильно изменяется. Если в "традиционном" случае в растворе превалирует ион РЬС142" ([РЬС142"] = 4.806-10~4). то в варианте "точного" расчета его концентрация снижается почти на три порядка ÜPbCl42"J = 8.748-Ю"7), а основными формами миграции свинца становятся РЬС12° ([РЬС12°]=2.181-Ю"5) и РЬС13" ([РЬС13"]=1.204-Ю"5). Другим весьма важным следствием смены формы миграции синца является уменьшение растворимости PbS в случае "точного" расчета при 400°С ( mPbS = 3.477-Ю"5 ) по сравнению с "традиционным" ( mPbS = 5.428-Ю"4 Крис. 4).

При снижении температуры по мере "выкипания" С0г величина диэлектрической проницаемости смеси приближается к значению для чистой Н20, и распределение РЬ по формам приближается к "традиционному" (рис. 4). Уже при 300°С вновь превалирующей формой становится РЬС142". а растворимость галенита в "точном" варианте расчёта (mPbS = 2. 328-10"4) довольно близка к "традиционной" (гаРЬ5 = 2.693-Ю"4). В результате кривая растворимости PbS имеет резко выраженный максимум при Т ~ 275°С (рис 46).

Полученный результат позволяет иначе представить себе механизм

рН

181] -2 -4 -6 -8

а)

рН

6 5 4

3

ГС

200 300 400

Дв С1

-2

-4 -6 -8

РЬ2+

6 5 4 3

Т°С

200 300 400

оэ

Рис.4. Растворимость галенита шРЬБ [моль/кг], равновесные моляльные концентрации свинца и его хло-ридных комплексов и рН раствора для принятой в работе модели эволюции флюида. Варианты: а - "традиционного" и б - "точного" расчетов.

рудоотложения из вскипающего флюида. На первой стадии при остывании (400 * 275°С) растворимость РЬБ растет, так что флюид при своем движении будет все более насыщаться свинцом, растворяя рассеянные в породе зерна галенита, начиная с 275°С растворимость РЬЗ начинает снижаться и. следовательно, из богатого РЬ флюида начнет выпадать РЬБ. образуя локализованные рудные тела. Налицо эффект переотложения и концентрирования. Интересно, что полученные из расчёта температуры интенсивного рудоотложения (200-275°С) достаточно близки к общепринятым значениям, полученным из анализа газово-жидких включений. Существенно также и то. что наблюдаемый эффект переотложения будет иметь место только при эволюции богатых СОг флюидов, снижение диэлектрической проницаемости в которых заметно. Тем самым становится понятным, почему такие флюиды являются "рудообразующими".

В заключение сформулируем защищаемые положения.

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

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

3. Получена новая термодинамическая информация для компонентов систем Аб'Ц-Б^-О-Н и БЬ11 '-Б1 '-О-Н. Полученные значения термодинамических параметров компонентов позволяют достичь адекватного описания поведения систем для всего набора экспериментальных данных. Ис-

пользование уравнение состояния HKF дает основание для надежной экстраполяции их термодинамических свойств до 500°С и 3 кбар.

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

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

Содержание диссертации опубликовано в следующих работах:

1. Акинфиев H.H. Алгоритм расчета гетерогенных равновесий на микро-ЭВМ. Геохимия, 1986. N 6. с. 882-890.

2. Акинфиев H.H., Демина О.В.. Никоноров А. П.. Панкратов A.B.. Покровский В.А.. Умрихин В.А. Пакет программ для расчета гетерогенных равновесий на микро-ЭВМ. Тезисы докладов II Всесоюзного симпозиума Термодинамика в геологии. 1988. г. Миасс. т. 1, с. 14-15.

3. Акинфиев H.H., Никоноров A.n., Панкратов A.B.. Умрихин В.А. Термодинамический анализ равновесий в системе As(III)-S(II)-0-H. Тезисы докладов И Всесоюзного симпозиума Термобинашка е геологии. 1988, г. Миасс, т. 1. с. 16.

4. Акинфиев Н.Н.. Никоноров А. П.. Умрихин В.А. Пакет программ для расчета фазовых равновесий на микро-ЭВМ. Тезисы докладов II Всесоюзного совещания Физико-химическое моделирование в геохиши и петрологии на ЭВМ. 1988, г. Иркутск, ч. 1. с. 137.

5. Акинфиев Н.Н.. Никоноров А. П.. Осмоловский И. С. Методические указания по выполнению термодинамических расчетов компонентного состава природных вод. МГРИ, 1988. 28 с.

6. Акинфиев Н.Н., Никоноров А.П.. Осмоловский И. С. Методические указания к лабораторным работам по теме Термодинамические расчёты в растворах. МГРИ, 1989. 22 с.

7. Акинфиев Н.Н., Никоноров А.П.. Осмоловский И.С. Оценка кинетики процесса подземного выщелачивания на основе физико-химического моделирования с использованием ЭВМ. Доклады координационного Совета по подземному выщелачиванию. Изд-во ВНИИКТ. 1989, фонды МГРИ. N Ф-1370.

8. Акинфиев Н.Н.. Никоноров А.П. Модель кинетики растворения минералов. Геохимия. 1989, N 6, с. 897-906.

9. Акинфиев Н. Н.. Келебеев А. С., Никоноров А. П., Умрихин В.А. Моделирование процессов кристаллизации и растворения аурипигмента на ЭВМ. Тезисы докладов IV Всесоюзной конференции по массовой кристаллизации и кристаллизации методом разделения смесей. Иваново. 1990. с. 72.

10. Акинфиев Н.Н. Дбмина О.В.. Никоноров А.П.. Панкратов А.В.. Умрихин В. А. Методические указания по термодинамическим расчетам равновесий в растворах на ЭВМ. МГРИ. 1991, 28 с.

11. Акинфиев Н.Н. Демина О.В.. Епифанова С.С.. Никоноров А.П.. Панкратов А. В., Умрихин В.А. Методические указания для расчёта полей устойчивости минералов в координатах р-Т. МГРИ, 1991. 26 с.

12. Келебеев А. С., Акинфиев Н. Н., Панкратов А. В.. Никоноров А. П.. Евжанов X., Умрихин В.А. Оптимизация выделения стронция из глубинных вод Туркмении. Тезисы докладов конференции МГРИ Новые материалы в области наук о Земле. Москва. 1991, с. 96-97.

13. Акинфиев Н.Н., Зотов А. В.. Никоноров А.П. Термодинамический анализ равновесий в системе As(III)-S(II)-0-H. Геохимия, 1992. N 5. с. 721-734.

14. Aklnflev N.N., Shikina N. D.. ZotovA. V. Thermodynamic properties of aqua species of Sb(III)-0-H-S(II) system at 25 -500°C and 1-1500 bar. Abstracts of the Second Int. Symp. Thermodynamics of Natural Processes, 1992, Novosibirsk, p. 5.

15. Aklnflev N.N. Calculation of thermodynamic properties of ion pairs at high state parameters. Abstracts of the Second Int. Symp. Thermodynamics of Natural Processes. 1992. Novosibirsk, p. 76.

16. Акинфиев Н.Н. Расчёт термодинамических свойств ионных пар при высоких параметрах состояния. Геохимия. 1992. N 11, с. 1416-1425.

17. Акинфиев Н.Н.. Зотов А. В.. Шикина Н.Д. Экспериментальные исследования и согласование термодинамических данных в системе Sb(III) - S(II) - 0 - Н. Геохимия. 1993. N 12. с. 1709-1723.

18. Акинфиев Н.Н. Программа БАЛАНС для расчета равновесий б мулъти-спстемах. Тезисы докладов конференции МГРИ поеие достижения е науках о Земле. Москва. 1993, с. 70.

19. Акинфиев Н. Н.. Епифанова С. С. Келебеев А. С.. Панкратов А. В.. Швец В. М. О гидрогеохимическом способе защиты от парникового эффекта. Геоэкология (Инженерная геология. Гидрогеология. Геокрино-логия). 1993. N 6, с. 18-26.

20. Акинфиев Н. Н.. Умрихин В. А. Алгоритм и программа SOL для расчета равновесий в природных и техногенных водах. Сборник рефератов докладов международной научной конференции Геофизика и современный тр. Москва. 1993, с. 243-244.

21. Aklnflev N. N.. Eplfanova S.S.. Kelebeev A. S., Mamrenko A. V.. Pankratov A. V.. Shvets V.M. Hydrogeochemlcal protection against greenhouse effect. Proceedings of Industrial and Agricultural Impacts on the Hydrologic Environment. Alexandria. 1993. p. 3.

22. Акинфиев H.H. Расчет термодинамических свойств водных частиц при повышенных температурах и давлениях. Парциальные мольные объемы электролитов. Геохимия. 1994, N 5. с. 681-690.

23. Акинфиев Н. Н.. Умрихин В. А. Алгоритм и программа SOL для расчета равновесий в природных и техногенных водах. Геоинформатика. 1994, N 2(2), с. 60-66.

24. Aklnflev N.N. Algorithms for ab initio calculation of ion association In supercritical H20 fluids. Abstracts of the 13th IUPAC Conference on Chemical Thermodynamics. 1994, Clermont-Ferrand. France, p. 133-134.

25. Швец В. M.. Акинфиев Н. Н.. Осмоловский И. С., Щербаков А. С. Теоретическое обоснование усовершенствования локального мониторинга качества подземных вод. В кн. Геология, ч. II/ Ред. кол.: Тихонов А. Н.. Садовничий В. А. и др. - М.: Из-во МГУ. 1994, с. 47-52.

26. Akinflev N.N. Model for computing mineral formation from boiling fluid: account of dielectric constant. Abstracts of the 16th General Meeting of the International Mineralogical Association, 1994, Pisa, Italy, p. 5-6.

27. Акинфиев H.H. Модель для расчета рудоотложения из вскипающего флюида: учёт диэлектрической проницаемости. Геохимия, 1994, N 10, С. 1445-1456.

28. Акинфиев Н.Н. Модель для расчета ab initio ионной ассоциации в сверхкритических водных флюидах: 1:1 электролиты. Геохимия. 1994, в печати.

Работа проводилась при финансовой поддержке проекта "Геомодель"

програшы "Университеты России".