Моделирование диффузии

[Elvis777]

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

 

Вопрос как узнать, как будет меняться концентрация в каждой точке с течением времени.

 

Компьютерная модель будет выглядеть следующим образом. Пространство сосуда будет разбито на кубики. Для простоты представим кубик 3x3x3(всего 27 кубиков) как в кубике Рубика. В самом начале, каждый из кубиков заполнен водой.

Затем, в один прекрасный момент один из заполненных водой кубиков (срединный прилегающий к одной из сторон кубика Рубика) оказывается заполненным медом. То есть в 26 кубиках концентрация меда 0%, а в одном 100%. И концентрация начинает выравниваться.

 

(Шаг первый)

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

 

(Шаг второй)

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

 

 

Вопросы:

Правильно ли я составил модель? (шаг один и шаг два)

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

Где можно прочитать описание диффузии с той точки зрения с которой мне это нужно.

 

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

 

Мне кажется должно получится что-то вроде этого.

http://www.google.ru...0%B2%D0%B0%D0%BD%D0%B8%D0%B5&um=1&hl=ru&newwindow=1&sa=N&biw=1280&bih=643&tbm=isch&tbnid=56retuUYrulz7M:&imgrefurl=http://freelance.ru/users/Dimitris/%3Fwork%3D199854&docid=1lrU_AxypVmqjM&itg=1&imgurl=http://freelance.ru/img/portfolio/pics/00/03/0C/199854.jpg&w=760&h=516&ei=kPrAT7O6DsTvsga2y-SjCg&zoom=1&iact=hc&vpx=790&vpy=161&dur=127&hovh=185&hovw=273&tx=94&ty=117&sig=114318136980987483479&page=1&tbnh=135&tbnw=199&start=0&ndsp=19&ved=1t:429,r:5,s:0,i:80

 

http://www.google.ru...0%B2%D0%B0%D0%BD%D0%B8%D0%B5&um=1&hl=ru&newwindow=1&sa=N&biw=1280&bih=643&tbm=isch&tbnid=chDFj74u8en5WM:&imgrefurl=http://www.philosophy.nsc.ru/journals/philscience/16_03/00_egorov.htm&docid=tvtD1QGrpWxQ6M&imgurl=http://www.philosophy.nsc.ru/journals/philscience/16_03/00_egorov.files/image002.jpg&w=623&h=426&ei=kPrAT7O6DsTvsga2y-SjCg&zoom=1&iact=hc&vpx=677&vpy=316&dur=531&hovh=186&hovw=272&tx=120&ty=102&sig=114318136980987483479&page=1&tbnh=129&tbnw=188&start=0&ndsp=19&ved=1t:429,r:10,s:0,i:91

Изменено 26 Мая, 2012 в 18:31 пользователем Elvis777

[arkansas]

Правильно ли я составил модель? (шаг один и шаг два)

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

Где можно прочитать описание диффузии с той точки зрения с которой мне это нужно.

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

Лучше не вставлять такие длинные ссылки, они все равно не отображаются полностью, и посмотреть удается с трудом.

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

О законах диффузии можно прочитать хотя бы в википедии http://en.wikipedia.org/wiki/Fick's_laws_of_diffusion

Ваш наглядный пример с медом так и остался для меня непонятен. Что такое нейромедиаторы, тут все знают. И с кинетикой изменения концентраций знакомы. Чай не кулинары все-таки.

Позволю себе спросить: кто поставил такую задачу? И еще. Вы с математическими пакетами (типа Maple, Matlab, Mathematica) знакомы?

Изменено 26 Мая, 2012 в 23:15 пользователем arkansas

[Wergilius]

Как заметил arkansas, в случае простой диффузии, ставится начально-краевая задача для уравнения Фика, которая решается одним из численных методов (с помощью конечно-разностных схем).

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

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

Изменено 27 Мая, 2012 в 15:46 пользователем Wergilius

[Elvis777]

Посмотрел я на уравнение Фика точнее я его и ранее видел но я не понимаю как его можно применить в той ситуации, которую нужно решить мне.

в вики http://uk.wikipedia....%96%D0%BA%D0%B0 написано следующее:

 

В системе с градиентом концентрации вещества dC/dx в направлении х диффузионный поток J определяется первым законом Фика

 

J = -DdC/dx,

где D — коэффициент диффузии (знак «-» указывает на направление потока от больших концентраций к меньшим).

 

В случае градиента концентрации не только в направлении х, нужно использовать общую формулу:

J = -D ∇ μ

где μ — химический потенциал.

 

2) В системе с градиентом концентрации вещества dC/dx в направлении х скорость изменения концентрации вещества в данной точке, обусловлена

 

диффузией, определяется вторым законом Фика:

 

dC/dt=Dd^2/dx^2,

где t — время.

 

Вроде химию учил немного, но ничего из этого не понимаю.

Главные вопросы: Как применить формулы? и Как составить алгоритм пересчета концентраций во всех кубиках? (о второй проблеме позже)

 

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

 

Представим одномерную последовательность кубиков (я говорю постоянно о кубиках по тому что я буду это программировать на основе массивов (одно двух и трехмерных), а кубики и массивы это одно и тоже, почти).

 

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

 

Итак, у нас есть одномерная последовательность кубиков. В левом боковом кубике (кубик №1) концентрация растворяемого вещества 100%, а в остальных, правее лежащих кубиках, концентрация растворителя 100%.

 

Первый шаг (цикла, итерации)

Теперь растворяемое вещество из кубика №1 попадает в кубик №2, а растворитель из кубика №2 в кубик №1.

 

Второй шаг(цикла, итерации) тут будет вопрос о алгоритме пересчета концентраций.

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

 

Но вопрос того как правильно просчитывать это второй вопрос. Первый и самый главный Как производить пересчет концентраций?

 

Упростим модель до двух кубиков. До сосуда, в котором есть перегородка и одна часть сосуда заполнена растворяемым веществом а другая растворителем на 100%. И в один прекрасный момент перегородка исчезает. И концентрация в разные моменты времени будет меняться. Вопрос: какой она будет через 1 через 2, 3, 4, 5 ... и т.д. секунд. Ведь сразу после исчезновения перегородки разница концентраций будет больше и значит, уравнивание концентраций будет происходить быстрее, чем в минуты перед самым завершением. То есть на первом шаге(как я понимаю) большее кол-во вещества перейдет из одного полупространства в другое чем на последнем шаге.

Как мне формуле второго закона Фика это вычислить(J = -DdC/dx, ).

И как это сделать с помощью следующей формулы J = -D ∇ μ

где μ — химический потенциал.

 

Не могли бы Вы привести пример для двух пространств (для первых 3-4 шагов)

 

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

 

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

 

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

 

То что для точного решения нужно сделать кубики очень крошечными это понятно(насколько точна модель будет зависеть от компьютера, тем более, что точность модели регулировать можно изменением пары параметров). Но я совсем понял фразу «концентрации за конечное время и на конечном расстоянии не выравниваются ни в каких из кубиков» Можно растолковать? Или это обратная сторона бесконечности?

Изменено 3 Июня, 2012 в 21:18 пользователем Elvis777

[Elvis777]

Вот так это должно выглядеть в графике.

[arkansas]

То что для точного решения нужно сделать кубики очень крошечными это понятно(насколько точна модель будет зависеть от компьютера, тем более, что точность модели регулировать можно изменением пары параметров). Но я совсем понял фразу «концентрации за конечное время и на конечном расстоянии не выравниваются ни в каких из кубиков» Можно растолковать? Или это обратная сторона бесконечности?

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

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

Я быстро решил одномерную задачу, используя стандартные функции Матлаба для численного решения диф. уравнений в частных производных и 2-й закон Фика. Начальным условием была произвольная протяженность начальной концентрации по координате. В качестве граничных условий было принято отсутствие диффузии на краях, т.е. условие Нейманна - равенство градиентов нулю. Преимуществом решения 2-го закона Фика является и то, что краевые условия можно максимально приблизить к реальным - например, учесть обратный захват или метаболизм нейромедиаторов, если кинетика известна из литературы. Так что, если целью не ставится изобретение велосипеда, считать диффузию методом кубиков (конечных разностей) имеет смысл лишь для двух- и трехмерных задач. Хотя я уверен, что в библиотеках современных инженерно-физических пакетов (типа Comsol) имеются средства расчета диффузии для набора разных геометрий и начально-краевых условий. Но если есть сугубое желание попрактиковаться в программировании или совсем плохо с соответствующим софтом, тогда другое дело...

Но вопрос того как правильно просчитывать это второй вопрос. Первый и самый главный Как производить пересчет концентраций?

Допустим, имеется трехмерный массив с начальной концентрацией C₀(x,y,z), в котором каждый элемент соответствует элементарной пространственной ячейке (кубику). Задан интервал времени (шаг) Δt и длина ребра кубика Δl.

Шаг первый. Расчет изменения концентрации за счет диффузии происходит для каждого из кубиков в соседние через три, или четыре, или пять, или шесть граней, в зависимости от пространственного расположения ячейки. Количество продиффундировавшего вещества через грань равно -DSΔtΣ_i(ΔC_i,i±1/Δi), где i=x,y,z. , S- поверхность грани, D - коэффициент диффузии. Поскольку мы зафиксировали, что Δi=Δl, неизвестной переменной будет только ΔC_i,i±1 - разность концентрации с прилегающими ячейками. Ее легко найти простым сдвигом исходного массива концентрации вправо и влево для каждого из трех направлений и вычитанием этого массива из полученных. Легко видеть, что изменившаяся после первого шага концентрация в кубике будет C₀(x,y,z)+DΔtΣ_i(ΔC_i,i±1)/(Δl)², где Σ_i(ΔC_i,i±1) - сумма разностей концентрации со всеми соседними кубиками

Второй шаг повторяет все операции первого с изменившейся начальной концентрацией. И т.д. до желаемого результата.

Только у меня есть серьезные сомнения, что эта задача в действительности моделирует какой-либо реальный эффект в синапсе, и вот по какой причине. Влияние адвекции/конвекции, наверное, можно не учитывать, принимая, что эффлюзия содержимого синаптического пузырька идет с малой скоростью. Но как быть с мембранными потенциалами, которые могут достигать несколько десятков мВ? Я несилен в нейрофизиологии, но мне кажется, что вакуоль не открывается при деполяризации мембраны. Принимая длину синаптической щели 20 нм, получаем электрическое поле в синапсе 100-1000 кВ/м. При рН около 7 очень многие эндо- и экзогенные нейромедиаторы, агонисты и антагонисты (за редкими исключениями нейтральных каннабиноидов, сальвинаринов) имеют положительный заряд. В таком поле будет развиваться стремительный электрофорез, который напрочь перебьет любую диффузию. ИМХО нервный импульс - настолько тонкая вещь, что вряд ли можно ожидать его регулирование только таким неспецифическим процессом, как диффузия. Гораздо вероятнее самосогласование через поляризацию/деполяризацию/реполяризацию и электромиграцию.