четверг, 16 декабря 2010 г.

О плагинах и линковке

В какой-то момент развития Indigo — кроссплатформенной C++ библиотеки для решения задач химической информатики — мы осознали, что любая нормальная библиотека такого рода должна позволять third-party расширения. Сказано-сделано. Indigo обрела поддержку плагинов. Обёртки на Java, Python и C# также расширяемы, при этом плагин на соответствующем языке является обёрткой бинарного модуля-плагина на C++. Всё это работает на Windows, Linux и Mac OS X. Дальше речь пойдёт о проблемах, с которыми мы столкнулись при линковке бинарных модулей, и о том как эти проблемы решить.

Постановка задачи

Есть библиотека на C++, назовём её «Core». В ней есть какое-то количество классов и некоторый интерфейс на C, который используется для «оборачивания» библиотеки в Java/Python/C# модули.

Есть ещё одна библиотека на C++, назовём её «Plugin». В ней тоже есть классы и C-интерфейс. Plugin имеет доступ к C++ классам Core и её глобальным функциям, тогда как Core не в курсе о существовании Plugin.

Core линкуется в динамическую библиотеку, назовём её libcore. Plugin тоже линкуется в динамическую библиотеку, назовём её libplugin. Эта библиотека не содержит в себе libcore.

Задача 1: libcore должна загружаться в адресное пространство программы на Java/Python/C# при инициализации соответствующей обёртки Core. Тип операционной системы (Windows/Linux/Mac OS X) и «битность» (32/64) должны определяться автоматически, т.е. нужно выбрать из коллекции сборок libcore для всех платформ правильную, и загрузить именно её. Сразу оговоримся, что обёртки поставляются в комплекте со всеми вариантами библиотеки, чтобы таким образом не нарушать «переносимость» и не причинять головную боль разработчикам, которые будут использовать обёртку.

Задача 2: libplugin должен загружаться в адресное пространство программы на Java/Python/C# при инициализации соответствующей обёртки Plugin. Само собой, и в этом случае ОС и битность определяются автоматически. Кроме того, libplugin должна «увидеть» уже загруженную libcore и подцепить из неё нужные C++ классы и глобальные функции.

Решение задачи 1

В Linux и gcc помогут ключи -m32 и -m64, с которыми собирается соответственно 32-битный и 64-битный код. Для кросс-компиляции, т.е. сборки под платформу, не соответствующую той, на которой выполняется компиляция, надо установить пакет gcc-multilib.

В Windows и Visual Studio дело делается установкой платформы Win32 или x64 в свойствах проекта. Чтобы платформа 64 появилась на 32-битных версиях VS, надо при инсталляции VS поставить галочку напротив «x64 Compilers and Tools», ну или доустановить потом отдельно.

На Mac OS X «битность» не является проблемой благодаря технологии универсальных бинарников. Не иначе как для компенсации этого удобства, бинарники, собранные для версии Mac OS X 10.6, не запускаются на 10.5, так что имеет смысл собирать для 10.5 и 10.6 отдельно. Или собирать только для 10.5 и рассчитывать, что они подойдут для 10.6, но мы это не пробовали.

Теперь пара слов о том, как определить платформу на этапе выполнения программы на Java/Python/C# и загрузить нужный модуль libcore.

Java

Узнать тип операционной системы можно с помощью System.getProperty("os.name"). «Битность», соответственно, через System.getProperty("os.arch"). Версия Mac OS X находится в System.getProperty("os.version").

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

Python

Константа os.name равна "nt", если дело происходит на Windows, и "posix", если это Linux или Mac OS X. Последняя отличается тем, что platform.mac_ver() возвращает определённую структуру данных, из которой как раз можно узнать версию системы. Битность можно узнать с помощью platform.architecture()

Загрузку бинарного модуля в Python, как и дальнейшую работу с ним, проще всего осуществлять с помощью ctypes.

C#

Операционная система в данном случае однозначно Windows. (Нет, до Mono мы ещё не добрались.) Битность можно определить, например, проверкой IntPtr.Size. Для вызова C-функций из C# используется DllImport("core.dll"). При вызове любой функции, объявленной с таким атрибутом, CLR попытается загрузить core.dll из системных директорий в память, если библиотека с таким именем ещё не загружена. Поскольку наша core.dll лежит не в системной директории, и вообще существует в двух вариантах (32 и 64 бита), то мы должны загрузить её до того, как будет вызвана любая функция из Core. Для того, чтобы загрузить core.dll из произольной директории, нужно воспользоваться системным вызовом LoadLibrary, который, в свою очередь, доступен через DllImport("kernel32").

Решение задачи 2

Функции, а также методы C++ классов, которые то же самое что функции, в терминологии линковщика называются символами. Символы, которые библиотека «показывает» наружу, называются экспортируемыми символами. Символы, которые библиотека берёт из других библиотек, называются импортируемыми символами. Надо сделать так, чтобы символы libcore, которые используются в libplugin, экспортировались из libcore и импортировались в libplugin. В каждой ОС динамический линковщик работает по-своему, и рецепты соответственно для каждой ОС свои.

Linux

С экспортом и импортом принципиальных сложностей нет. линковщик (ld), который собирает объектные файлы в .so-модуль, все неопределённые (undefined) в библиотеке глобальные символы считает импортируемыми, а все определённые символы считает экспортируемыми (есть возможность экспортировать только избранные символы, но нам это не нужно). Динамический линковщик (который не ld, а часть операционной системы) при загрузке libcore запишет все экспортированные из неё символы, а при загрузке libplugin обнаружит в ней символы, ждущие импорта, и импортирует их в libplugin из libcore.

Есть ещё одна особенность. Символы могут быть помечены импортируемыми «откуда угодно», или импортируемыми из конкретной библиотеки.

Импорт символов из конкретной библиотеки

ld помечает символ импортируемым из конкретной библиотеки если при линковке указать .so-файл, из которого этот самый символ экспортируется. То есть, в libplugin.so будет чётко прописано, что такой-то символ импортируется из библиотеки libcore.so и только из неё. Кроме того, динамический линковщик будет пытаться при загрузке libplugin.so загрузить libcore.so из директорий, записанных в:

  1. Переменной окружения LD_LIBRARY_PATH. При этом, учитывается только значение этой переменной на этапе запуска программы. Фокусы с подменой LD_LIBRARY_PATH по ходу пьесы не работают. Для программ с setuid/setgid LD_LIBRARY_PATH и вовсе игнорируется.
  2. Кэше /etc/ld.so.cache
  3. /lib
  4. /usr/lib
(man ld-linux)

Если libcore.so не будет найдена ни в одной из указанных директорий, то загрузка libplugin.so не пройдёт успешно. Нетрудно понять, что для наших целей такой подход не годится, т.к. libcore мы распространяем в двух вариантах (32 и 64 бита) и обязательно вместе с самой программой, чтобы разработчики на Java и Python не терпели неудобств с непереносимостью своих программ из-за бинарных модулей.

Импорт символов из любой библиотеки

Если в ld при линковке libplugin.so не передавать libcore.so, то он пометит отсутвующие символы как импортируемые, но не укажет откуда именно. Динамический линковщик затем при загрузке libplugin.so не станет пытаться загрузить libcore.so, а попытается найти отсутствующий символ среди всех загруженных на данный момент библиотек (+в самой программе). Конечно, libcore.so будет среди них, т.к. мы инициализировали Core до того, как начали инициализацию Plugin. Всё очень хорошо, но есть ещё одна деталь.

На самом деле, ld перебирает не все загруженные на данный момент библиотеки, а только те из них, которые были загружены с флагом RTLD_GLOBAL (man dlopen). Те библиотеки, которые загружены с флагом RTLD_LOCAL, подходят только для импорта символов конкретно из них (см. предыдущий заголовок). Виртуальная машина Java и её System.load(), что бы вы думали, конечно загружает все библиотеки с RTLD_LOCAL, без вариантов. Но есть обходной путь! Он появился в glibc 2.2 (2000 г.) не иначе как специально для Java: это флаг RTLD_NOLOAD. После вызова System.load("/path/to/libcore.so") можно вызвать (уже не из Java, а из C):

dlopen("/path/to/libcore.so", RTLD_NOLOAD | RTLD_GLOBAL);
и символы из libcore.so, ранее «закреплённые» за libcore.so, станут доступны для импорта по схеме «откуда угодно». Динамический линковщик отдаст их в нужный момент в libplugin. Дополнительная изюминка заключается в том, что вызов dlopen для libcore.so с RTLD_NOLOAD можно оформить в самой libcore.so, в какой-нибудь функции инициализации. Будет работать.

Что касается загрузки .so-файлов в Python, то с ним всё гораздо проще, т.к. ctypes поддерживает произвольные флаги для загрузки библиотек, в т.ч. и RTLD_GLOBAL:

lib = CDLL("/path/to/libcore.so", mode=RTLD_GLOBAL);

Mac OS X

Правила линковки на Mac OS X такие же, как и в линуксе, есть только небольшие различия в терминологии и опциях линковщика. Схема с привязкой символов к библиотеке, описанная в предыдущем разделе, здесь имеет название: two-level namespace (man ld). Альтернативная схема, которая нам и нужна, называется flat namespace. На этапе компиляции libplugin надо выбрать, по какой схеме испортировать в неё символы. Для нас это означает, что надо передать в ld ключ "-flat_namespace".

Использование flat namespace, по замыслу разработчиков, не отменяет необходимости задания в командной строке ld всех зависимых .dylib-файлов (т.е. в нашем случае libcore.dylib). Это странно, но оправдано для более сложных случаев с косвенными зависимостями (indirect dynamic libraries), которые к нашей задаче не имеют отношения. Можно, тем не менее, попросить ld закрыть глаза на неразрешённые зависимости, передав ему ключ "-undefined suppress". В этом случае динамический линковщик будет разрешать их в рантайме, как в Linux.

Трюк с RTLD_NOLOAD | RTLD_GLOBAL в Mac OS X тоже актуален.

Windows

В Windows нет динамического линковщика.

Его там не может быть в принципе из-за отсутствия поддержки position-independent code, но от этого не легче. В Windows есть только загрузчик динамических библиотек (loader), который, мягко говоря, не совсем в курсе про динамическое связывание. Вся работа по динамическому связыванию делается самой загружаемой DLL на этапе инициализации. Компилятор и линковщик (link.exe) должны позаботиться о том, чтобы DLL сделала эту работу правильно. Программист, в свою очередь, должен позаботиться о том, чтобы компилятор и линковщик правильно поняли свою задачу.

Экспортируемые и импортируемые символы

По указанным выше причинам, в DLL не может быть неопределённых (undefined) символов. Никто не проверит на этапе загрузки DLL, какие символы в ней «defined», а какие «undefined». Такого вопроса просто не стоит. Все символы должны быть определены. В том числе и символы, которые импортируются из другой DLL.

Этот парадокс разрешается следующим образом: при компиляции модулей DLL, в которую импортируются символы (в нашем случае plugin.dll) компилятор на месте импортируемых символов создаёт функцию-"прослойку", которая

  1. Загружает в память DLL, в которой находится нужная функция (LoadLibrary("core.dll")). Большая удача, нет мороки с путями. Есди загрузить core.dll из нужной директори заранее (LoadLibrary("\path\to\core.dll")), то plugin.dll не станет искать core.dll в системных директориях, а просто «подхватит» уже загруженную копию.
  2. Получив указатель (handle) на загруженную DLL, ищет в ней нужную функцию (GetProcAddress) и вызывает её с теми параметрами, которые были переданы в неё саму, т.е. в прослойку.

link.exe не догадается сделать такую прослойку для всех C++ функций, которые не определёны в DLL (хотя догадается сделать это для C-функций, слабое утешение). Перед объявлением каждой из C++ функций, которую вы хотите импортировать в plugin.dll, надо писать __declspeс(dllimport).

Более того, link.exe не догадается экспортировать из DLL символы, которые в ней определены. Перед каждой функцией, которую вы хотите экспортировать из core.dll, надо писать __declspec(dllexport). Это касается и C, и C++ функций. При линковке любой DLL помимо собственно DLL возникает маленький файл с расширением .lib — он и содержит вышеуказанные «прослойки» для всех экспортированных функций.

Получается, что одни и те же C++ функции должны объявляться в Core как __declspec(dllexport), а в Plugin — как __declspec(dllimport). Объявляются эти функции в заголовочных файлах, которые одинаково включаются в Core и в Plugin. Самое время воспользоваться препроцессором:

#ifdef _WIN32
#ifdef PLUGIN
#define DLLEXPORT __declspec(dllimport)
#else
#define DLLEXPORT __declspec(dllexport)
#endif
#else
#define DLLEXPORT
#endif
и перед всеми функциями Core, которые нужны в Plugin, писать макрос DLLEXPORT. При сборке plugin.dll надо указать компилятору макрос PLUGIN.

Экспортирование C++ классов

Физически, C++ класс в скомпилированных модулях — это набор его методов, т.е. функций. Экспортирование класса означает экспортирование всех его методов, кроме приватных. Экспорт и импорт классов осуществляется тем же образом, что экспорт и импорт функций. Можно использовать тот же макрос DLLEXPORT, что и для функций.

Написав DLLEXPORT класса, у которого есть суперкласс, или неприватные поля — объекты каких-либо классов, или неприватные методы, возвращающие объекты каких-либо классов, то при компиляции вы заметите следующий варнинг:

warning C4251: *** : class *** needs to have dll-interface to be used by clients of class ***
Это, в принципе, правильное предупреждение. Если Plugin импортирует класс «CoreClass», определённый в Core, то он конечно будет использовать его публичные (public) методы и публичные поля. Или отнаследуется и будет использовать защищённые (protected) методы и поля. Все классы, возникающие в процессе работы с CoreClass, должны быть доступны через импорт так же, как и сам CoreClass. И всех их надо тоже пометить как DLLEXPORT, о том и речь в варнинге C4251. Иначе Plugin ждёт ошибка линковки.

Если бы в проектах Core и Plugin не было бы классов-шаблонов, то на этом бы наше повествование благополучно закончилось. Но шаблоны у нас есть. Нет, сами классы, которые надо экспортировать, не являются шаблонами, но среди их методов есть такие, которые принимают или возвращают объекты шаблонных классов. Вот на тему этих шаблонных классов и возникает C4251.

Можно, конечно, попытаться писать DLLEXPORT при объявлении шаблонов. В простых случаях это даже будет работать. Но на самом деле это лишено всякой логики. Шаблон — это не класс, он никогда не экспортируется. Экспортируется конкретный класс, который создаётся в тот момент, когда компилятор встречает инстанцированный шаблон. Выйдет так, что все классы-экземпляры «экспортируемого» шаблона будут экспортироваться из Core и импортироваться в Plugin. Это хорошо до тех пор, пока Plugin не воспользуется одним из шаблонов Core с параметром, которого нет в Core. (То есть, из Core попросту не экспортируется данный экземпляр шаблона). И даже тогда всё может быть в порядке; но рано или поздно окажется, что компилятор не в состоянии понять, что не надо импортировать данный экземпляр шаблона из Core, а надо его создать. В нашем случае это произошло, когда один «экспортированный» шаблон с неэкспортированным экземпляром инстанцировался в другом, тоже с неэкспортированным экземпляром. Дело закончилось ошибкой вида
error LNK2001: unresolved external symbol "__declspec(dllimport) ..."

Есть ещё один способ обойти C4251, под названием «explicit template instantiation», т.е. явное экспортирование экземпляра шаблона. Он подробно разбирается в этой статье, в применении к STL, с которой кстати приписать DLLEXPORT к объявлению шаблона нет возможности. Код получается немыслимо громоздким, а результаты — неутешительными. В статье отмечен тот факт, что экспортирование шаблонных классов в случае линковки нескольких (независимых друг от друга) DLL, содержащих одни и те же шаблонные классы, приведёт к ошибкам линковки вида

LNK1169: one or more multiply defined symbols found
Да, link.exe, в отличие своего коллеги ld, не умеет определять и отбрасывать дубликаты функций в разных библиотеках.

Короче говоря, экспортировать шаблоны нет смысла и вообще нельзя. Есть два выхода из ситуации с C4251:

  1. Игнорировать, благо это варнинг, а не ошибка. core.dll и plugin.dll просто будут иметь по пачке одинаковых методов для шаблонных классов. Конфликта при загрузке не будет, поскольку классы не экспортированы. Ошибки линковки тоже не будет, если все не-шаблонные классы имеют в своём заголовке DLLEXPORT.
  2. Не экспортировать свои классы, а экспортировать только их методы, с помощью того же DLLEXPORT в объявлении каждого метода. Как ни странно, варнинг при этом исчезает. Странно, потому что опасность-то остаётся. Опасность заключается в том, что core.dll и plugin.dll всё равно будут иметь независимые реализации одних и тех же шаблонных классов; и если окажется так, что plugin.dll была собрана с иной версией этих классов, бинарно несовместимой с той, с которой была собрана core.dll, то они не смогут вместе работать с этими классами. В статье по ссылке выше сказано, что именно это может произойти с классами STL, когда одна библиотека собрана под VS7, а другая под VS8.

Прочие проблемы с Windows

Ещё пара неприятностей, по числу которых Windows и так лидирует с большим отрывом:

  1. Чтобы собрать plugin.dll, надо поставить в проекте Plugin зависимость на проект Core. Однако, если появилась необходимость завести в проектах Plugin и Core конфигурации сборки «Static», и собирать их в plugin-static.lib и core-static.lib, то зависимость останется, и её не убрать. core-static.lib будет вкомпиливаться в plugin-static.lib, вместе со своими глобальными переменными. Проект, который зависит от обеих библиотке Core и Plugin, при сборке обретёт две копии Core, и работать скорее всего не будет из-за наличия двух копий глобальных данных. Придётся создавать отдельные проекты CoreStatic и PluginStatic. Даже тогда, не получится оставить Core и Plugin пустыми и просто поставить им зависимости соответственно от CoreStatic и PluginStatic: при линковке Core и Plugin не будут экспортироваться символы. Придётся дублировать в Core/Plugin и CoreStatic/PluginStatic одинаковые наборы файлов.
  2. В компиляторе VS есть опция /MT, что означает, что рантайм (CRT) будет «вкомпилен» в библиотеку или программу. Это очень удобно, с учётом того что даже в самые современные версии Windows эта CRT не включена! В Windows есть только очень старая версия, а VS линкует с новой, которая доступна для Windows либо в составе VS, либо в виде отдельного пакета. Так вот, опция /MT в случае двух связанных между собой DLL приведёт к ошибке во времени выполнения. Она неизбежно случится, когда например память выделена в одной DLL, а освобождается (или даже копируется) другой DLL. Остаётся использовать опцию /MD и обеспечивать затем присутствие «Redistributable Package» на машине пользователя.

вторник, 13 июля 2010 г.

Навигация в мире органических соединений

Сколько органических соединений вы знаете? А сколько вы знаете лекарств? Каждое лекарство, не считая тех, что производятся из растений, преставляет собой комплект из действующего вещества и оболочки/растворителя, в которой/ом оно проходит свой путь до усваивания организмом пациента. Действующее вещество в лекарственном препарате — это одно конкретное органическое соединение1. Количество известных органических соединений, которые можно добыть или синтезировать, превышает 30 миллионов, а количество лекарств на рынке — всего несколько тысяч. Создание любого нового лекарства занимает от 10 до 15 лет и является очень дорогостоящим. Расходы на программное обеспечение составляют в этой индустрии (как и почти в любой другой) весьма скромную долю.

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

Роль компьютеров в этом процессе за последние два десятилетия стремительно возросла, и вот почему. Производство лекарственного средства — комплексная задача, в которой есть место пробам и ошибкам.

Представим, что усилиями биологов в организме выявлен «нездоровый» белок, вызывающий болезнь или болезненные ощущения. Дело за малым — найти вещество, которое разрушит или заблокирует белок, не причинив вреда организму в процессе. Затруднение состоит в том, что на эту роль может годиться одна молекула из 30 миллионов. Или ещё не открытая молекула. Современные технологии массового синтеза (т.н. комбинаторная, или сочетательная химия) и массовых биохимических опытов (HTS), позволяют за короткие сроки получать сотни тысяч новых молекул и гигабайты экспериментальных данных.

Опишем карьеру лекарственного вещества от конца к началу. До того, как попасть на прилавки аптек, лекарство должно пройти клинические испытания (clinical trials) на пациентах, под присмотром врачей. Это представляет определённый риск, поэтому до пациентов доходят лишь немногие вещества, прошедшие доклинические тесты (preclinical testing) на животных. Их тоже берегут, поэтому доклиническим испытаниям предшествуют массовые опыты на отдельных живых клетках (in vitro, лат. «в стекле», т.е. в пробирке). Но и в пробирки не бросаются все молекулы подряд. Люди должны выбирать нужные (перспективные для лечения данной болезни) вещества и отбрасывать заведомо неподходящие, и без компьютера им этого не сделать2. На самом деле, и с ним не очень удобно. Эффективная система навигации по химическим содинениям пока ещё не создана; и о перспективах создания таковой сейчас пойдёт речь3.

Поиск химических соединений в базах данных

Состояние систем поиска химических соединений в наши дни, увы, примерно соответствует состоянию поисковых систем во всемирной паутине в 90-е годы4. Да, именно так. Примитивные алгоритмы поиска (и поиск ведётся далеко не по всем имеющимся источникам); весьма вялая поддержка языковой грамматики; никакого ранжирования результатов.

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

Например, амфетамин является подструктурой мезокарба, и оба препарата подавляют реакцию обратного захвата дофамина (dopamin reuptake) в мозге, что ведёт к повышению активности. Но мезокарб, в отличие от амфетамина, не вызывает тахикардии и повышения артериального давления.

Амфетамин и мезокарб

Вообще говоря, нередко случается так, что добавление или удаление небольшого фрагмента идёт молекуле (точнее, пациентам) на пользу. Чем добавлять и удалять всевозможные фрагменты вручную, проще запустить поиск по базе данных и найти все молекулы, содержащие данную как подструктуру или все молекулы, содержащиеся в данной. Соответствующие виды поиска называются «подструктурный поиск» (substructure search) и «надструктурный поиск» (superstructure search).

Подструктурный поиск в сервисе Bingo с амфетамином в качестве запроса.

Более общий критерий структурного сходства молекул основан на количестве различных фрагментов, которые присутствуют одновременно в обеих молекулах. Поиск молекул по такому критерию назывется поиском по сходству (similarity search).

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

Найденные в базе данных молекулы не очень интересны сами по себе; они интересны в контексте. В каких медицинских и химических статьях упомянуто данное соединение? Есть ли на него патент? Представлено ли оно в коммерческих каталогах? Известны ли его свойства, такие как растворимость, кислотность, токсичность, внутренняя абсорбция и другие? Известна ли химическая реакция его синтеза? Доступны ли исходные компоненты этой химической реакции? На эти вопросы и на многие другие должна давать ответ система поиска.

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

  • PubChem — база данных из 27 миллионов соединений с богатыми возможностями для поиска: по номеру, по названию, по структурной формуле, по подструктуре и по сходству. Химические свойства также можно задавать в качестве дополнительных критериев поиска (например, ограничиться только молекулами, молекулярная масса которых не превышает 120). Базу постоянно пополняют более 80 организаций.
  • ChemSpider содержит 25 миллионов соединений и имеет важное отличие от PubChem: добавлять молекулы и обновлять информацию о них здесь могут не только избранные организации, но и простые пользователи. Вместе с последними, список источников ChemSpider составляет почти 300 пунктов. Поиск соединений в ChemSpider не имеет такого количества опций, как в PubChem; в частности, отсутствует поиск по сходству.
  • eMolecules — компиляция из 7 миллионов соединений, собранных в 150 коммерческих каталогах. Возможности поиска минимальны; никакой информации о соединениях, кроме ссылок на каталоги, сайт не показывает. Это скорее платформа для продавцов химических веществ, нежели поисковая система для исследователей.

Поиск химических соединений по научным публикациям.

Пионерами поиска в научных работах по химии были создатели химической реферативной службы (Chemical Abstracts Service, CAS), существующей с 1907 года. В этой службе ведётся учёт всех известных химических соединений. Тысячи людей в течение десятков лет вручную составляют библиографические справки и заполняют базу данных SciFinder, отдельного продукта CAS для поиска публикаций. Аналогичная база данных, поддерживаемая издательством Elsevier, называется «Crossfire Beilstein». Сервисы PubChem и ChemSpider также выдают пользователю вместе с каждой найденной молекулой список публикаций, к которым данная молекула может иметь отношение; но возможности для поиска собственно публикаций в этих сервисах не очень развиты.

SciFinder

Как же, наконец, отправить на отдых тысячи «индексаторов», от рассвета до заката читающих статьи и заполняющих библиографические базы? Эта задача несколько сложнее, чем найти слово «парацетамол» по текстам статей. Во-первых, само вещество может иметь несколько названий (пример альтернативного названия парацетамола — «N-(4-гидроксифенил)ацетамид»). Во-вторых, лекарство, содержащие это вещество, может упоминаться под разными торговыми марками (в данном случае «Панадол», «Эффералган» и десяток других). В третьих, вещество может быть не написано, а нарисовано в статье. В растровом виде (в старых отсканированных статьях), или в векторном (начиная с 90-х годов). Программы по автоматическому распознаванию рисунков с молекулами сегодня находятся в плачевном состоянии. Вот наиболее известные из них:

  • CLiDE канадской фирмы SimBioSys
  • OSRA — проект с открытыми исходниками нашего соотечественника, работающего в США
  • ChemoCR— проект Марка Циммермана из института Фраунгофера в Германии
CLiDE — наиболее развитая программа из перечисленных, но она нередко ошибается, требует вмешательства человека и «не знает» многих особенностей молекул. OSRA активно развивается, но обладает на данный момент худшим качеством распознавания. ChemoCR, похоже, находится в перманентной закрытой разработке: эту программу никто никогда не видел, тем не менее доступно немалое количество публикаций по алгоритмам, используемых в ней. Указанные программы ещё менее пригодны к распознаванию более сложных химических объектов, как-то: химических реакций, таблиц с заместителями. Комбинированный семантический анализ текста и рисунков (например, «молекула на рис. 10a имеет показатель LD50 равный 5.6 г/кг для взрослых крыс») вообще нигде не реализован.

Планирование синтеза

Схема синтеза химического соединения представляет собой цепочку химических реакций.

Многоступенчатая реакция синтеза парацетамола

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

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

Перегруппировка Брука
Перегруппировка Кляйзена
Благодаря этому обстоятельству, значительно увеличивается разнообразие синтезируемых веществ, и в то же время усложняется поиск. Вкупе с «многоступенчатостью», планирование синтеза становится сложной задачей. В некоторой степени эта задача решена в сервисе Reaxys, который, увы, доступен только по подписке.

Reaxys

Заключение

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

В настоящий момент наиболее перспективная система, которая может в будущем стать «Гуглом и Википедией» химиков — это ChemSpider. Её создатели явно принимают в расчёт ключевые элементы успеха глобальных сервисов: кросс-платформенность (работает через броузер и даже с мобильных устройств), доступность для каждого, богатые возможности, «дружба» с многими другими сервисами (включая PubChem), право пользователей публиковать свой контент. ChemSpider имеет недостатки, как общие для индустрии, так и свои собственные, но движется в правильном направлении.

ChemSpider

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

Благодарности

Автор признателен Н. Велецкому и Д. Лушникову за ценные замечания по тексту статьи.

Сноски

1 В редких лекарствах их два, например в бисептоле.

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

3 На прочих стадиях роль вычислительных машин не менее важна; однако возникающие там задачи лежат за рамками данной статьи.

4 (до появления Google). Тогдашних «королей» информационного поиска (Altavista, Lycos, Rambler) мало кто помнит; в том числе оттого, что они были практически бесполезны.

5 Взаимодействие молекулы с белком тоже можно моделировать на компьютере. Этому посвящена отдельная область вычислительной химии под английским названием «docking».

Примеры неудовлетворительной работы OSRA

(Картинки выложены по просьбе Игоря в комментариях). На первых двух и на последней OSRA не выдаёт никакого результата, на 3-й, 4-й и 5-й выдаёт неверный результат.

воскресенье, 30 мая 2010 г.

Symbolic Integration

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

Классы функций одной переменной

Чем шире класс функции, тем сложнее брать интеграл. Например, интегрирование полиномов одной переменной
p(x) = an xn + ... + a1 x + a0
не составляет никакого труда. Уже некоторую сложность представляет интегрирование рациональной функции
r(x) = an xn + ... + a1 x + a0

bm xm + ... + b1 x + b0
.
Если степень числителя n меньше степени знаменателя m, и коэффициенты лежат в поле вещественных чисел, то интегрирование прозводится следующим образом: знаменатель раскладывается на множители вида (x - xi)ki и (x2 + bj x + cj)lj, где каждый xi — корень кратности ki, а bj и cj задают пару сопряжённых комплексных корней кратности lj. Затем методом неопределённых коэффициентов дробь r(x) раскладывается на сумму конечного количества дробей вида
Ai

(x-xi)ki
и Bj x+Cj

(x2 +bjx + cj)lj
Для каждого из трёх видов слагаемых (трёх, т.к. дробь с Bjx + Cj разбивается на отдельные слагаемые с Bjx и Cj) интеграл находится с помощью формул из справочника. В этих формулах встречаются логарифмы, радикалы, модули и арктангенсы, поэтому интеграл рациональной функции в общем случае не является рациональной функцией.

Если степень числителя n больше или равна степени знаменателя m, то мы можем разделить числитель на знаменатель и получить сумму полинома и «правильной» дроби. Если нас интересуют функции с комплексными коэффициентами, то всё делается точно так же, как и с вещественными, с тем отличием, что квадратичные члены из разложения знаменателя будут упразднены.

Кажется, что интегрирование рациональных функций — уже чисто техническая задача? Это не совсем так. Если степень многочлена больше 4, то его корни в общем случае не записываются в радикалах. Это доказал Абель в 1826 г. Например, вещественный корень уравнения «x5 + x + a = 0» невозможно записать, пользуясь арифметическими операциями и операцией извлечения корня.1

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

Алгебраическая функция y от переменной x задаётся неявным образом, через полином от x и y:
pn(x) yn + ... + p1(x) y + p0(x) = 0
Для примера преобразуем указанное выше уравнение 5-й степени в алгебраическую функцию y от переменной x:
y5+y+x=0
Алгебраические числа и функции тоже не являют собой предел выразительности математической записи. Хорошо известны трансцендентные числа (e, π) и трансцендентные функции (экспонента, логарифм, тригонометрические функции), которые невозможно представить через полиномиальные уравнения.

К классу элементарных функций относятся алгебраические функции, экспонента, логарифм и любые их конечные композиции. Подразумевается оперирование над полем комплексных чисел, что означает, что элементарны и все тригонометрические функции, т.к. они выражаются через комбинации eix и e-ix. Обратные тригонометрические функции, соответственно, выражаются через логарифмы. Показательная функция ax выражается через экспоненту и логарифм: ax = ex ln(a).

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

Интегрирование в теории

С рациональными функциями, несмотря на указанные выше неприятности, принципиальных трудностей нет: интеграл всегда существует в виде элементарной функции, и есть процедура, которая его гарантированно найдёт.

С алгебраическими функциями дело обстоит сложнее.

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

Это даёт следующую общую запись для интегрируемых алгебраических функций:
f(x) = v0'(x) + m
Σ
i=1
ci vi'(x)

vi(x)
(1)

Кроме того, Лаплас утверждал, что интеграл функции можно записать так, чтобы он содержал лишь те алгебраические величины, которые есть в исходном подынтегральном выражении.

Гипотезу Лапласа доказал Лиувилль в 1833 г., а ещё через два года он опубликовал доказательство аналогичного факта для некоторого обобщения класса алгебраических функций, т.н. элементарных функций по Лиувиллю. Они имеют форму алгебраических функций, но в роли переменных в них не только x, но ещё конечный набор логарифмов и экспонент от алгебраических функий от x:

y(x) = f(x, z1(x), ... zr(x)) — алгебраическая функция от x, z1(x), ... zr(x)
zi(x) = ln(qi(x)) или exp(qi(x)), где qi — алгебраическая функция от x.

Лиувилль показал, что если элементарная (по Лиувиллю) функция имеет элементарную первообразную, то эта первообразная представляется следующим образом:

f(x, z1(x), ... zr(x)) dx = y0(x, z1(x), ... zr(x)) + m
Σ
i=1
Ai ln(yi(x, z1(x), ... zr(x)))
Выражаясь менее формальным языком, первообразная функции y — это алгебраическая функция от того, что уже имеется в y, плюс конечное количество логарифмов (с коэффициентами) алгебраических функций того, что уже имеется в y. То есть, принцип Лапласа остаётся в силе.

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

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

Итак, дифференциальная теория Галуа. Ключевое понятие в ней — это т.н. элементарное расширение поля функций2, которое получается последовательным добавлением в него логарифмов, экспонент и алгебраических функций:
  • ln(f(x)), где f ∈ K
  • exp(f(x)), где f ∈ K
  • алгебраической функции, полученной уравнением вида P(f1(x), ... fn(x))=0, где f1, ... fn ∈ K.
(за K считается поле, расширенное предыдущим добавлением). Нетрудно видеть, что любая элементарная функция принадлежит какому-либо элементарному расширению поля рациональных функций.

Ключевая теорема дифференциальной теории Галуа является обобщением теоремы Лиувилля и говорит о том, что если функция f из поля F имеет первообразную g в поле G — элементарном расширении F, то f записывается в форме (1), где vi ∈ F, а g соответственно содержит функцию из F и логарифмы функций из F с весами.

Если строить поле F «специально» под функцию f, т.е. подбирать минимальное расширение поля рациональных функций, содержащее F, то станет очевидно, что принцип Лапласа снова не нарушается: первообразная функции f содержит то, что уже имеется в f, а также логарифмы того, что уже имеется в f.

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

Интегрирование на практике

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

Первая реализация алгоритма Риша была выполнена Джоэлом Мозесом в рамках знаменитого «Project MAC» в MIT в 1971 г. Программа под названием SIN интегрировала чисто трансцендентные функции.

Дэвенпорт в 1981 г., основываясь на работе Риша и некоторых глубоких результатах дифференциальной алгебры и комплексного анализа, разработал алгоритм интегрирования чисто алгебраических функций и реализовал его в известной среде символьных вычислений REDUCE-2.

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

Алгоритм Дэвенпорта ещё и обладал большой вычислительной сложностью, что его автор признавал и выражал надежду, что алгоритм можно упростить. Так и произошло: Барри Трагер из MIT в 1984 г. внёс серьёзные улучшения в алгоритм Дэвенпорта. Обновлённый алгоритм обладал гораздо большей практической ценностью и был реализован в математических программах Axiom и Maple.

Решающий шаг к практическому решению вопроса сделал Мануэль Бронштейн в 1990 г, обобщив алгоритм Трагера на произвольные элементарные функции.

В 1998 г. Бронштейн написал монографию по символьному интегрированию, в которой понятным языком изложил лучшие достижения в этой области, начиная с теоремы Лиувилля и заканчивая собственными результатами.

Бронштейн был одним из ведущих разработчиков Axiom. Вместе с Трагером он реализовал в Axiom свой алгоритм, но лишь частично. В 2005 году Бронштейн скончался, спустя два года после смерти основателя Axiom Ричарда Дженкса.

При своей неполноте, имплементация Бронштейна — наиболее полная из имеющихся в программах символьных вычислений. Автор отмечал, что при возникновении нереализованного случая Axiom выдаёт ошибку, в отличие от процедур интегрирования Maple и Mathematica. То есть, если Axiom обнаруживает, что интеграл от данной функции не является элементарной функцией, то это правда. Если Axiom выдаёт ошибку, то интеграл может быть или не быть элементарной функцией.

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


x4 + 10 x2 - 96 x - 71
равна
- 1

8
ln((x6+15 x4-80 x3+27 x2-528 x+781)

x4+10 x2-96 x-71
- x8 - 20 x6 + 128 x5 - 54 x4 + 1408 x3 - 3124 x2 - 10001),
однако в таком виде её находит только Axiom. Прочие пакеты — Mathematica, Maple, Matlab, Maxima, Reduce — в лучшем случае выдают громоздкий ответ с участием эллиптических интегралов, т.е. неспособны выдать ответ в виде элементарной функции.3

Первообразную функции earcsin(x), равную ½(√ 1-x2 +x)earcsin(x), способны найти Mathematica, Axiom и Reduce.

Неопределённый интеграл ∫√arctan(x) dx невозможно найти (или определить его отсутствие) ни в одном из перечисленных пакетов; при этом Axiom выводит наиболее осмысленное сообщение: implementation incomplete (constant residues).

Нерешённые вопросы

Полная имплементация алгоритма Риша-Трагера-Бронштейна ждёт своего автора. Однако, и это ещё не всё: в теории символьного интегрирования тоже остались открытые вопросы.

Пришло время сказать, что алгоритм Риша (даже для чисто трансцендентных функций) не является «честным» алгоритмом.

Дело в том, что в процессе работы алгоритм Риша выстраивает «башню» элементарных расширений поля рациональных функций, с тем чтобы «верхушка» башни содержала подынтегральное выражение. При этом не всегда понятно, готова ли уже «башня» или нужно добавить в неё ещё один логарифм, экспоненту, или алгебраическую функцию. Например, если в подынтегральном выражении встретилась экспонента ex, она добавляется в поле (происходит элементарное расширение поля функцией ex). Если при дальнейшем анализе в подынтегральном выражении встречается ex2, то эту функцию также следует добавить в поле, вырастив таким образом «башню» до двух этажей. Но если после этого попадётся функция ex+x2, то алгоритм должен определить, что она принадлежит сформированному на данный момент элементарному расширению, и не наращивать «башню». Да, в описанном случае это несложно сделать; чуть сложнее определить, например, на языке элементарных функций, что arctan(x) = -arctan(1/x). Но в общем случае эта задача очень сложна4, и нынешние реализации алгоритма Риша используют лишь эвристики для её решения, которые по определению ненадёжны.

Как указывает Дэвенпорт и следом за ним Бронштейн, при реализации алгоритма можно «сделать вид», что нам известен т.н. базис трансцендентности для поля констант, то есть например считать e и √ 2 алгебраически независимыми. Это значит, что мы принимаем на веру, что нет никакого полинома от e и √ 2 , который был бы равен 0. Такое решение нужно принимать о каждой новой константе, появляющейся в ходе работы алгоритма, которую эвристика не позволяет связать алгебраическим выражением с ранее учтёнными константами. Если полученный базис трансцендентности окажется неверным, это может отразиться на корректности выданного ответа в том случае, если ответ гласит, что «интеграл не является элементарной функцией» (да, это касается и Axiom, о которой выше было сказано противоположное). Если же интеграл найден, его можно считать правильным.

Известным частным случаем задачи об алгебраической независимости является т.н. Constant Problem: определение того, тождественно ли нулю некоторое математическое выражение. Задача тривиально решается для многочленов, но для более широкого класса функций, как показал Даниэль Ричардсон в 1968 г., эта задача неразрешима.

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

Если есть функция f, получающаяся с помощью сложения, умножения и композиции из синуса, экспоненты и модуля, а коэффициенты f используют рациональные числа, π и ln(2), то задача проверки утверждения «f тождественно равна нулю» алгоритмически неразрешима.

Функция модуля существенна в доказательстве Ричардсона, и есть основания полагать (да и сам автор теоремы надеется), что без неё задача Constant Problem всё таки разрешима, что даёт надежду на создание абсолютно честного алгоритма символьного интегрирования.

И последнее: класс элементарных функций кажется несколько надуманным, не так ли? Например, функция ошибки
erf(x) = 2



π

x

0
e-t2 dt
не является элементарной, как доказал Лиувилль. Но чем она хуже экспоненты и логарифма, которые мы относим к элементарным функциям? На компьютере она может быть вычислена с произвольной точностью, и даже, как экспонента и логарифм, обладает элементарной производной. Эту функцию вполне можно включить в «надэлементарные» расширения функций, так же как и эллиптические интегралы, ∫xxdx и многие другие. Тема интегрируемости в «надэлементарных» функциях была затронута ещё автором SIN в 1970-х годах, и ныне представляет обширное поле для исследований.

Благодарности

Автор признателен А. В. Демьянову, П. А. Мозоляко и Н. Н. Васильеву за вычитку черновиков этой статьи и высказанные замечания.

Литература

J. Liouville. Premier memoire sur la determination des integrales dont la valeur est algebrique. Journal de l’Ecole Polytechnique, 14:124–148, 1833

J. Liouville. Second memoire sur la détermination des integrales dont la valeur est algebrique. Journal de l’Ecole Polytechnique, 14:149–193, 1833

J. Liouville. Mèmoire sur l’intégration d’une classe de fonctions transcendantes. J. Reine Angew. Math. Bd. 13, p. 93-118. (1835)

G. H. Hardy. The Integration of Functions of a Single Variable, 1916.

R. H. Risch (1970). The solution of the Problem of Integration in Finite Terms. Bulletin of the American Mathematical Society 76 (3): 605–608.

Joel Moses. Symbolic integration the stormy decade. Proceedings of the second ACM symposium on Symbolic and algebraic manipulation (1971) 427-440.

Дж. Дэвенпорт. Интегрирование алгебраических функций, Москва «Мир» 1985.

B. Trager. Integration of algebraic functions. PhD Thesis, MIT, 1984

M. Bronstein. Integration of elementary functions. Journal of Symbolic Computation, Volume 9, Issue 2 (1990). 117-173

M. Bronstein. Symbolic Integration Tutorial. ISSAC'98, Rostock.

Richardson, D. Some Unsolvable Problems Involving Elementary Functions of a Real Variable. J. Symbolic Logic 33, 514-520, 1968.

Сноски:

1Кстати, это число называется корнем Бринга Br(a), и с его помощью можно записать все комплексные корни более общего уравнения «x5 + px + q = 0».
2Формально говоря, дифференциальная теория Галуа оперирует дифференциальными полями характеристики 0, но мы в эти термины погружаться не станем и будем считать по-прежнему, что имеем дело с различными классами функций.
3Проверялись следующие версии программ: Axiom 2010, Mathematica 7, Maple 13, Matlab R2007, Maxima 5.21, Reduce 2009.
4Например, до сих пор никто не знает, является ли e + π рациональным числом (хотя известно, что e и π по отдельности — трансцентентные числа). Известно только равенство ei π = -1 и бесконечная дробь Рамануджана, связывающая e и π.

воскресенье, 21 февраля 2010 г.

Mark Tarver и Qi

Qi (произносится «чи») — функциональный язык программирования, с проверщиком типов полным по Тьюрингу, как и в Haskell и C++. Автор (Марк Тарвер, профессор в отставке) при этом утверждает, что его творение затмевает все существующие разработки и останется непревзойдённым навсегда, по причине использования наиболее гибкой нотации для описания типов. К тому же, оно работает поверх Лиспа и разрешает использовать функции на Лиспе.

О практическом применении языка речи не идёт; на вопрос, а что такого есть в мире написанного на Qi, автор гордо отвечает, что сам Qi на себе и написан, не считая вспомогательных библиотек к нему же. Скорее, Qi — исследовательский инструмент для разработки новых систем типов; ну и для саморазвития, конечно.

По части последнего Марк вообще увлечённая личность: в августе прошлого года помпезно попрощался с Qi и вообще с программированием, и уехал в Индию заниматься йогой. В ноябре посетовал, что-де в ашраме туго с интернетом, а в декабре (уже не с такой помпой) оповестил сообщество, что «phase of travels is over», вернулся в Англию и продолжает развивать Qi.

Как тут не пожелать доктору Тарверу всяческих успехов.