photon190573: (Default)
По просьбам населения попробую разложить по полочкам, что такое мультиреференсные методы и как они соотносятся друг с другом.

В HF/DFT волновая функция представлена в виде 1 детерминанта Слейтера. 1 детерминант (или 1 спин-адаптированная комбинация детерминантов) описывает 1 электронную конфигурацию. Соответственно, если состояние (неважно, основное или возбужденное) не описывается 1 конфигурацией (например, частично заполненная оболочка или ситуации, когда состояния примерно одной энергии описываются разными эл. конфигурациями), такое состояние требует описания в виде лин. комбинации спин-адаптированных волновых функций, соответствующих разным конфигурациям.

Например, состояние 2p^1 -- это 3 конфигурации: 2p_x^1, 2p_y^1, 2p_z^1. Если не учесть их все с равными весами, волновая функция состояния будет перекошенной (например, будет только p_x в сферически симметричной задаче, что сходу неверно).

Соответственно, любое описание состояния через лин. комбинацию конфигураций, является многоконфигурационным (да, в CIS или TDDFT возб. состояния тоже являются многоконфигурационными, даже если доминирует одна). И не так уж важно, получено ли это описание вариационным методом (любые методы ССП) или по теории возмущений, или через кластерное разложение.

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

Но Full CI, если это не очень маленькая система -- это непозволительно дорого. Количество конфигураций растет как факториал (быстрее, чем экспонента). Хочется как-то облегчить себе жизнь. А давайте ограничимся только 1-2-кратными возбуждениями по отношению к референсной ("возбуждения" при этом могут иметь ту же энергию что и референс) -- это CISD (КВ1+2). Про CIS тут даже и говорить как-то неприлично... Другой вариант: а давайте учтем все возможные конфигурации, но на ограниченном наборе орбиталей (например, в пределах 1 оболочки). Это КВ в полном активном пространстве (CAS-CI). Полное -- потому что учтены все возможные для данных активных орбиталей расположения электронов. А активные -- это те самые орбитали, к-рые мы для этих целей выбрали. Ну и можно совсем уж урезать осетра и сделать ограниченное КВ в акт. пространстве. Правда, сейчас так никто не делает за бессмысленностью.

Итак, если речь идет про КВ (CI), то орбитали не трогают, только оптимизируют коэф. КВ, и варианты бывают Full-CI, CAS-CI, CISD и CISD в акт. пространстве.

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

Что значит -- подправлять орбитали? Значит, что получив коэффициенты КВ и построив лин. комбинацию конфигураций, подставим ее в ур. Хартри-Фока и снова диагонализуем методом ССП, но с учетом того, что теперь у нас не 1 детерминант, а целая лин. комбинация. Вот тут и всплывает 1й нюанс: если мы, допустим, для нашей системы 2p^1 будем подправлять орбитали только для 1 состояния (к-рое у нас в ХФ получилось 2p_x^1), то оно и будет у нас получаться дальше. Это называется State-Specific (SS-MCSCF) Как было косое решение в симметричной задаче, так и осталось. Отсюда растут ноги у идеи усреднять по состояниям -- усредняется как энергия, так и матрица плотности, и индивидуальные конфигурации уже не играют такой роли (State-Averaged, SA-MCSCF). В случае вырожденных состояний этот прием отлично работает. Но и в случае квазивырождения, и даже вообще невырождения (когда между усредняемыми состояниями есть заметная энерг. щель) тоже очень неплохо работает: в одном расчете получаются сразу несколько состояний с одинаковой погрешностью. Если в SS-MCSCF каждое состояние получалось в отдельном расчете (только низшее состояние имеет смысл, а остальные хоть и получаются, но смысла особого не имеют), с риском свалиться в вариационный коллапс или с необходимостью вводить симметрийные ограничения, то в SA-MCSCF усредняемые состояния сразу ортогональны друг другу, и никакого вариационного коллапса. И симметрию можно отключить, особенно если и изначально молекула не была симметричной.

Итак, результат вариационного расчета называется референсом. Многоконфигурационный референс получается в результате MCSCF, в отличие от одноконфигурационного (из HF или DFT). MCSCF бывает без усреднения по состояниям (SS-MCSCF) или с усреднением (SA-MCSCF), бывает полным в акт. пространстве (CASSCF) или ограниченным (обозначений не встречала, хотя задать такой расчет в Гамессе или FF легко). Варианта CISD с оптимизацией орбиталей не встречала, хотя технически, наверняка, возможно. Но скорее всего, будет нереально дорого. Похоже что фишка многоконфигурационного референса как раз в том, чтобы получить решение в небольшом акт. пространстве, а потом его улучшать дальше.
Далее, имея на руках многоконфигурационный референс, можно (и нужно!) улучшать его еще дальше за счет взаимодействия состояний из акт. пространства с состояниями из т.н. внешнего пространства -- теми, где участвуют орбитали, не входящие в состав активных. Это все методы с приставкой MR: MRCI (обычно в варианте CISD) над референсом из акт. пространства, MRPT и MRCC -- соответственно, теория возмущений и Coupled Clusters на референсе из акт. пространства. Из них только MRCI -- вариационный.
Грабли там заботливо разложены на каждом шагу, но это уже отдельная история.
Размерность и прожорливость каждого метода -- это, наверное, не ко мне. Лимитировать могут разные вещи. Небольшая (в смысле размера базиса) задачка в огромном акт. пространстве порождает чудовищное кол-во конфигураций, отчего затыкается любой традиционный алгоритм (но специально заточенный под такие вещи DMRG, наверное, нет). Не особо большой SA-CASSCF(8e,8o) для задачи в огромном базисе тоже может стать кластеру поперек горла. Считали нормально CASSCF для N состояний, увеличили до N+1 -- и уперлись в нехватку памяти или безумное время счета.
photon190573: (Default)
Общение с пользователем катализовало меня на пост о целесообразности использования многоконфигурационных методов. В принципе, я об этом уже писала тут, но, видимо, стоит поподробнее. Я просто скопирую из лички свой ответ этому пользователю с некоторой редакцией.

Иногда в литературе попадаются странные вещи. Например, простое одноконфигурационное основное состояние исследуется методом CASSCF или ваще CASPT2. Смотришь на такое чудо мысли и думаешь: автор, зачем ты машину гонял казеную? (ну не на лаптопе же личном, в самом деле, CASPT гонять?) что ты поймать хотел? просто хотел продемонстрировать, что знаешь такое слово -- CASSCF? а поскольку взять акт. пространство приличного размера кишка оказалась тонка (ну или кластер слабоват :) ), то сделали CAS(2,2) -- 2 электрона на 2 орбиталях.
Пользователю, к-рый только начал осваивать многоконфигурационные методы, такое простительно. А автору статьи в журнале просто не попался въедливый рецензент и не заставил объяснять, зачем все-таки многоконфигурационные методы применили к одноконфигурационной задаче.
Итак, когда многоконфигурационные методы нужны, и когда -- нет.

Обычно к многоконфигурационным методам прибегают, когда либо основное состояние подозревается на неодноконфигурационность (открытые оболочки, где N электронов рассаживаются по M квазивырожденным орбиталям я-собьюсь-со-счета-сколькими способами), либо интересует описание основного и пачки возбужденных состояний на одном уровне. Либо исследуют ППЭ, где ожидаются участки с квазивырождением (т.е., фактически сразу придется исследовать целую пачку близколежащих ППЭ). В этих случаях используют усреднение по целевым состояниям.
Использовать многоконфигурационные методы только для основного состояния в случаях без квазивырождения -- стрелять из пушки по воробьям и с большой вероятностью мазать. Если интересует только основное состояние, то усреднение не делают -- а зачем? И тогда получается, что сначала мы в CASSCF слегка улучшим основное состояние добавкой некоторых конфигураций из акт. пространства, а оно обычно достаточно маленькое, и добавка даже десятка конфигураций с весом >=1% (реально там десятка и не наберется, 3-5 штук помимо самой доминирующей, а остальное -- на уровне пыли) погоды не делает, и вообще не факт, что возбуждения именно с этих активных орбиталей существенно улучшат описание корреляции в основном состоянии. И, не удовлетворившись этим, мы дополнительно улучшаем описание основного состояния во 2 порядке ТВ -- делаем расчет QDPT.
В таких случаях показано старое доброе MP2 -- относительно дешево, сердито, и корреляция должным образом учтена, и никакой лишней работы, все в одном флаконе.

Осмысленность совсем уж минималистских акт. пространств (типа (2,2) или (4,4)) тоже, имхо, в основном сомнительная.
Обычно такие акт. пространства берут, когда точно знают, что все квазивырождение обусловлено именно этими электронами именно на этих орбиталях. Например, в случае квазивырождения в открытой оболочке имеет смысл включать в акт. пространство только ее и ничего более. Тогда эл. состояния, возникающие в ней, будут описаны достаточно примитивно -- только в виде необходимых лин. комбинаций конфигураций из этого жалкого акт. пространства, но вся корреляция будет прекрасненько учтена во 2 порядке ТВ.
Другое дело -- когда исследуют возб. состояния. Там обычно априори ничего не ясно. Даже если CIS или TDDFT показывают, что в первом возбужденном состоянии доминирует HOMO->LUMO, какие еще возбуждения дадут вклад в его описание в CASSCF -- хз. Там обычно начинают с маленького акт. пространства (например, (2,2), т.е. HOMO+LUMO, но не обязательно, можно сразу больше брать) и понемножку расширяют. Причем добавлять орбитали по одной (сначала +1 к занятым, потом +1 к виртуальным и т.п.) -- плохая идея, т.к. такие акт. пространства обычно несбалансированы и либо сходимости не будет, либо сойдется к перекошенному решению. Как расширять акт. пространство -- вопрос проб и ошибок и приходит с опытом. Критерием служит заселенность нат. орбиталей (ее показывает Кемкрафт): близкая к 0 или 2 (~0.01 или 1.99, соответственно) заселенность означает, что данная орбиталь вклада в описание усредняемых состояний не дает, и ее спокойно можно выкинуть... иногда. Потому что если исследуется ППЭ, то орбиталь, не важная на одном участке, может оказаться ключевой на другом. Выкинем ее здесь -- не сможем гладко, без разрывов построить кривую "дотуда".

Так что если вдруг в литературе попался CASPT расчет с акт. пространством типа HOMO+LUMO, и при этом ни возб. состояния не интересуют, ни квазивырождения в основном состоянии не предвидится -- перед Вами типичное пускание пыли в глаза. Ну и сами так никогда не делайте :)
photon190573: (Default)
В прошлый раз я забыла написать про NSTATE в $XMCQDPT. Как я говорила, это размер эффективного гамильтониана.

Что такое эффективные операторы, можно почитать тут: http://www.chem.msu.ru/rus/teaching/zajzevskii/zaizev.pdf (лучше прочитать, чтобы понимать терминологию). Модельное пространство в MCQDPT натянуто на выбранные вектора КВ, полученные (обычно) в результате SA-CASSCF расчета (т.е., как раз те, что нас интересуют: они входили в wstate(1) в CASSCF и входят в wstate(1) и avecoe(1) в QDPT). Эти вектора -- только часть векторов CAS-CI, полное число к-рых равно числу конфигураций в CAS-CI. Если акт пространство небольшое, то полное кол-во конфигураций, соответствующих заданной мультиплетности, будет вполне обозримым. Так, в CAS(2,2) будет всего три синглетных конфигурации (ψA*ψA, ψB*ψB, (ψAα*ψBβ - ψBα*ψAβ)) и одна триплетная (ψA*ψB); в (CAS(4,4) синглетных конфигураций будет 20, а триплетных 15; а в CAS(6,7) септетных конфигураций будет 7. Конечно, в этом случае бессмысленно мудрить и брать размер модельного пространства меньше, чем полное число векторов CAS-CI. Но кол-во состояний в CAS-CI растет очень быстро, угнаться за ним в NSTATE невозможно. Поэтому для больших (>50 конфигураций) акт. пространств надо ограничивать размер эффективного гамильтониана каким-то разумным числом. Поскольку после всех преобразований (к-рые могут занять ооооочень долгое время) длительность дальнейшего QDPT расчета пропорциональна кол-ву состояний в NSTATE, то размер эффективного гамильтониана лимитируется либо таймлимитом Вашего кластера, либо разумным желанием не делать лишнюю работу.

Опыт показывает, что энергия первого возбужденного состояния (S1 в органическом красителе, без вырождения по симметрии) перестает заметно меняться при NSTATE=10. Если нужно больше состояний, то NSTATE надо увеличить пропорционально. Реально больше 40 состояний не берут, 20-30 -- вполне нормально.

Еще (по мотивам вопросов на Кемпорте) одно замечание. Если нам нужны состояния разной мультиплетности, и мы их изначально получили в детерминантном (cistep=aldet) CASSCF расчете с pures=.f. (собс-но, так и надо делать, если нам нужен общий набор орбиталей для состояний разной мультиплетности), в QDPT смешать разные мультиплетности не получится. Надо на одном и том же наборе $VEC делать отдельные QDPT расчеты для каждой мультиплетности. Для этого надо задать соответствующий MULT в $CONTRL и pures=.t. в $DET либо задать MULT в $XMCQDPT.

Итак, XMCQDPT расчет закончился. Ищем в выдаче строчку  *** XMC-QDPT2 ENERGIES ***. В табличке после нее будут вполне понятные вещи: энергии (CAS-CI на векторах из $VEC) исходных состояний и их исправленные во втором порядке ТВ энергии. Также нам могут пригодиться моменты переходов (RADIATIVE TRANSITION MOMENT). Для обычного спектра поглощения нужны будут RADIATIVE TRANSITION MOMENT BETWEEN STATE #   1 AND STATE # (все остальные), но в выдаче также присуствуют моменты переходов между всеми состояниями (уж не знаю, насколько они могут пригодиться). Еще могут понадобиться ZERO-ORDER PROPERTIES FOR STATE # (т.е., натуральные орбитали и распределение заряда в нулевом порядке ТВ для каждого из состояний), в частности, TOTAL MULLIKEN AND LOWDIN ATOMIC POPULATIONS (малликеновские заселенности).

Кое в чем нам поможет ChemCraft. Открыв им выдачу, мы можем увидеть непосредственно электронный спектр, где энергии переходов выдаются в нм, см-1 или эВ, а силы осциллятора вычислены из моментов переходов по известной формуле. Очень удобно. Еще ChemCraft покажет нам, переход с каких орбиталей на какие доминирует в каждом из возбужденных состояний. Tools->Orbitals->Render molecular orbitals -- и мы увидим натуральные (активные) орбитали с заселенностями в каждом состоянии. Заселенности, сильно отличающиеся от 0 или 2, как раз и указывают на орбитали, доминирующие в переходе. А распределение зарядов в каждом состоянии придется вытаскивать руками из TOTAL MULLIKEN AND LOWDIN ATOMIC POPULATIONS, поскольку ChemCraft нам тут не помощник.

ChemCraft показывает форму натуральных орбиталей (Zero order QDPT natural orbitals) для каждого из состояний в NSTATE. То, что для каждого из состояний форма одной и той же орбитали оказывается разной, -- нормально. Это в SA-CASSCF (из-за усреднения) орбитали были общими для всех состояний. Теперь мы каждое состояние уточнили индивидуально -- и вид орбиталей для каждого из них может стать своим собственным.

Напомню еще раз, что State 1 в CASSCF и QDPT -- это основное состояние данной (в MULT в $CONTRL) мультиплетности. Поэтому если нам нужны энергии, скажем, триплетов относительно основного синглетного состояния, то в выдаче для триплетов на энергии переходов, выдаваемые ChemCraft-ом (и отсчитанные от State 1 в данной выдаче, т.е. T1), смотреть не надо, а надо вычислить их ручками -- вычесть из них энергию S0 из выдачи для синглетов.

Некоторые тонкости.

Если посмотреть по строчке EIGENVECTORS OF THE EFFECTIVE HAMILTONIAN, то можно увидеть матрицу эффективного гамильтониана в базисе собственных векторов CAS-CI. Это то, во что превращается CAS-овский гамильтониан (в базисе этих векторов он диагональный) после включения возмущения. А возмущением здесь служат 1- и 2-кратные возбуждения во внешнее пространство. Если вид этого эффективного гамильтониана не сильно отличается от диагонального (на диагонали числа, близкие к 1, вне диагонали -- близкие к 0), то у нас вообще идеальная ситуация -- CAS-овский гамильтониан почти не придется уточнять по ТВ. Такого, конечно, не бывает. Поэтому очень хорошей считается ситуация, когда на диагонали стоят числа ~0.8-0.9, а вне диагонали числа ~0.1 и меньше. Но и так бывает редко. В обычном MCQDPT им. Накано большие внедиагональные элементы приводят к Большой Лаже. В XMCQDPT это не страшно, если перемешивание векторов происходит в пределах массивов состояний, задаваемых wstate(1) и avecoe(1) (т.е., в пределах модельного пространства). Если заметно примешиваются состояния, лежащие за пределами wstate(1) и avecoe(1), это повод получше посмотреть на исходный CASSCF и, скорее всего, что-то в нем серьезно скорректировать и пересчитать.

Если посмотреть по строчке EIGENVALUES OF THE NON-SYMMETRIC EFFECTIVE HAMILTONIAN, то там будет массив (вообще говоря, комплексных) чисел. Если мнимая часть (второй столбец) маленькая (отличная от 0 в 6 или более знаке), то все ОК. Если отличие от 0 появляется в знаке 5-м или меньше -- тоже повод пересмотреть исходный CASSCF.

Ну вот, вроде бы все. За кадром остался вопрос, для чего нам нужны одновременно состояния разной мультиплетности. Собс-но, понятно для чего: чтобы посчитать для них константы переходов за счет спин-орбитального взаимодействия. Но этот рассказ я оставлю на потом.
photon190573: (Default)

Итак, CASSCF сошелся. У нас есть чудненькие орбитальки, отлично согласующиеся с нашим представлением о системе. Ну или не очень согласующиеся, но не вызывающие резкого отторжения. Давайте посмотрим на энергии состояний, к-рые мы включали в усреднение. Нравятся? скорее всего, не очень. Иногда -- просто ни к черту. А что б вы хотели? это только начало!

Как я раньше уже говорила, ХФ не включает в себя электронную корреляцию вовсе, а МКССП включает только ее статическую часть. Динамическая корреляция (та, что связана с попытками электронов держаться друг от друга подальше) учитывается либо с помощью дополнительного КВ (уже над полученной многоконфигурационной волновой функцией), либо с помощью теории возмущений. На настоящий момент даже ограниченное (Singles, Singles+Doubles) КВ над CASSCF -- это дорого (хотя в отдельных случаях, на маленьких системах, его можно использовать в качестве бенчмарка), а теория возмущений 2 порядка -- вполне доступна.

Read more... )

XMCQDPT, реализованная в FireFly, представляет собой улучшение теории Накано, обладающее некоторыми полезными свойствами: (1) явная нетривиальная зависимость эффективного гамильтониана от размера модельного пространства в любом порядке ТВ; (2) сходимость энергий и др. свойств по отношению к расширению модельного пространства в каждом фиксированном порядке ТВ; (3) эффективный гамильтониан должен зависеть от всего подпространства, натянутого на на выбранные вектора CASCI, а не от конкретного выбора базиса в этом подпространстве; (4) полученные энергии должны быть однозначными, непрерывными и гладкими функциями от атомных координат и др. внешних параметров, с возможными исключениями в виде точек конических пересечений и др. областей случайного вырождения.


Как перейти от CASSCF к XMCQDPT? Во-первых, нужно оставить весь инпут от CASSCF, но в качестве группы $VEC взять OPTIMIZED ORBITALS от сошедшегося CASSCF расчета. Во-вторых, добавить ключ MPLEVL=2 в группу $CONTRL. В третьих, добавить следующие строчки:

 $XMCQDPT INORB=1 nstate=<кол-во состояний в эффективном гамильтониане>                                                 
  edshft=0.02 wstate(1)=<кол-во единичек, совпадающее с wstate в $DET или $GUGDM2>,-0 avecoe(1)=<кол-во единичек, совпадающее с wstate в $DET или $GUGDM2>,-0 $end
 $mcqfit $end

INORB=1 означает "не пересчитывать CASSCF, взять орбитали из $VEC и все настройки CASSCF из $MCSCF и $DET/$DRT, но проверить орбитали на ортогональность".
edshft -- величина т.н. сдвига энергетических знаменателей. Эта техника применяется для избавления от вторгающихся состояний. Вообще говоря, величина этого сдвига влияет на вычисляемые энергии, поэтому чем она меньше, тем лучше. Прелесть XMCQDPT в том, что для этого метода edshft=0.02 вполне достаточно, и трогать это значение не надо, тогда как в CASPT2, например, применяют значения ~0.3-0.4, и это как бы еще один способ "подрулить" энергии (т.е., еще один источник произвола). Присутствие группы $mcqfit $end включает некие полезные численные алгоритмы, значительно облегчающие жизнь.

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

В группу $SYSTEM добавим MWORDS=50 masmem=<намного больше>. Надо так распределить память, чтобы почти вся оперативная память досталась мастер-процессу (masmem), понемногу (50 мегаслов) получили бы остальные процессы, и сколько-то осталось бы системе. Понадобятся также строчки

$system mklnp=<кол-во процессоров на ноде> np=<кол-во процессоров на ноде> $end
$smp smppar=.t. $end

Остальные настройки -- такие же, как были в исходном CASSCF.

Продолжение следует. Анализируем результаты.

Profile

photon190573: (Default)
photon190573

April 2024

S M T W T F S
 123456
78910111213
14151617181920
2122232425 2627
282930    

Syndicate

RSS Atom

Style Credit

Expand Cut Tags

No cut tags
Page generated Jan. 17th, 2026 09:03 pm
Powered by Dreamwidth Studios