Уравнение, скриншот которого приведен на первом рисунке, является одним из вариантов модели ФитцХью-Нагумо. В этом уравнении переменная x отвечает трансмембранному потенциалу в биологической возбудимой ткани, а y — медленному току восстановления. Параметр ε задает временной масштаб между быстрой и медленной переменными; a является параметром возбуждения. В этой модели наблюдается состояние равновесия и возбудимый режим при a>1 и колебательный режим при a<1.
На рисунке 2 приведен типичный фазовый портрет нейрона ФитцХью-Нагумо в случае, когда параметр a задает колебательный режим.
На рисунке 3 приведен фазовый портрет для случая, когда уменьшен параметр временного масштаба между быстрой и медленной переменными (по сравнению с рисунком 2). Если на рисунке 2 вектора отрисовывались в течение 0.0055 единиц безрамерного времени, то на рисунке 3 — в течение 0.0015 единиц. При этом сами вектора выглядят горизонтальными, в отличие от случая ε=0.1. Это означает, что движения вдоль оси x значительно быстрее движений вдоль оси y, в случае малого параметр ε.
На рисунке 4 приведен пример фазового портрета для нейрона ФитцХью-Нагумо в возбудимом режиме, в котором наблюдаемым устойчивым состоянием является равновесие (черная точка). В целом направления движения траекторий в фазовом пространстве аналогично случаю, приведенному на рисунке 2, но отсутствует устойчивый предельный цикл, а появляется устойчивый узел, поэтому траектория останавливается в равновесии.
p. s. Чтобы сразу увидеть новый материал в моем блоге в своей ленте, подписывайтесь! Буду рад комментариям, вопросам, предложениям.
Программа для построения временной реализации нейрона ФитцХью-Нагумо
#include <iostream>
#include <fstream>
#include <cmath>
#include <iomanip>
#include <string>
#include <sstream>
#include <cstdlib>
#include <cstdio>
#include <ctime>
/*
/*The program written by Andrei Bukh, to compile it run:
* $ g++ -o FHN FHN.cpp -lm
* $ ./FHN 0.1 0.9 0.1 0.1 10000 8000 0.001
*/
int totalvars = 2;
int totalpars = 2;
int func (double *** vars, int marker, int srfr, int stto, double * times, double * pars) {
double x = vars[marker][srfr][0];
double y = vars[marker][srfr][1];
double eps = pars[0];
double a = pars[1];
vars[marker][stto][0] = (x*(1.0-x*x/3.0)-y)/eps;
vars[marker][stto][1] = x+a;
return 0;
}
int rk4step (double *** vars, double *times, double *pars, int marker, int maxdelay, double h) {
int pm = marker-1; if (pm<0) pm+=maxdelay;
for (int v=0;v<totalvars;v++)
vars[marker][5][v] = vars[marker][0][v] = vars[pm][5][v];
times[marker] = times[pm]; func (vars,marker,0,4,times,pars); /*K1*/
for (int v=0;v<totalvars;v++)
vars[marker][5][v] += h * vars[marker][4][v] / 6.0;
for (int v=0;v<totalvars;v++)
vars[marker][1][v] = vars[marker][0][v] + h * 0.5 * vars[marker][4][v];
times[marker] += 0.5*h; func (vars,marker,1,4,times,pars); /*K2*/
for (int v=0;v<totalvars;v++)
vars[marker][5][v] += h * vars[marker][4][v] / 3.0;
for (int v=0;v<totalvars;v++)
vars[marker][2][v] = vars[marker][0][v] + h * 0.5 * vars[marker][4][v];
func (vars,marker,2,4,times,pars); /*K3*/
for (int v=0;v<totalvars;v++)
vars[marker][5][v] += h * vars[marker][4][v] / 3.0;
for (int v=0;v<totalvars;v++)
vars[marker][3][v] = vars[marker][0][v] + h * vars[marker][4][v];
times[marker] += 0.5*h; func (vars,marker,3,4,times,pars); /*K4*/
for (int v=0;v<totalvars;v++)
vars[marker][5][v] += h * vars[marker][4][v] / 6.0;
return 0;
}
int singleparset (double *trvar, double *trpar, int maxdelay, double time, double trans, double timestep) {
double*** vars; double* times; double* pars;
vars = new double** [maxdelay];
for (int d=0;d<maxdelay;d++) vars[d] = new double* [6];
for (int d=0;d<maxdelay;d++)
for (int s=0;s<6;s++) vars[d][s] = new double [totalvars];
times = new double [maxdelay];
pars = new double [totalpars];
times[0] = 0.0;
for (int d=0;d<maxdelay;d++)
for (int s=0;s<6;s++)
for (int v=0;v<totalvars;v++) vars[d][s][v] = trvar[v];
for (int p=0;p<totalpars;p++) pars[p] = trpar[p];
std::ofstream ts_file("TS.dat");
std::ofstream ts_file2("TSvec.dat");
int marker = 0;
int t = 0;
ts_file2 << pars[0] << " " << pars[1] << " " << 0.0 << " " << vars[marker][0][0] << " " << vars[marker][0][1] << " ";
for (double ct=0.0;ct<time;ct+=timestep) {
t++;
marker++; if (marker>=maxdelay) marker = 0;
int prevmarker = marker-1; if (prevmarker<0) prevmarker = maxdelay-1;
rk4step (vars, times, pars, marker, maxdelay, timestep);
if (((t%(1000000))==0)) std::cout<<"ct="<<ct<<"\n";
if (ct>=trans) {
if (((t%(10000))==0)) std::cout<<"ct="<<ct<<"\n";
ts_file << pars[0] << " " << pars[1] << " " << ct << " " << vars[marker][0][0] << " " << vars[marker][0][1] << "\n";
}
}
ts_file2 << vars[marker][0][0] << " " << vars[marker][0][1] << "\n";
ts_file2.close();
ts_file.close();
for (int v=0;v<totalvars;v++) trvar[v] = vars[marker][0][v];
delete [] pars;
delete [] times;
for (int i=0;i<maxdelay;i++)
for (int j=0;j<6;j++) delete [] vars[i][j];
for (int i=0;i<maxdelay;i++) delete [] vars[i];
delete [] vars;
return 0;
}
float takefloatpar (int argc, char* argv[], int num) {
double res = 0.0;
if (argc > num) {
bool rightnumber = true;
try {
res = std::stod (argv[num]);
} catch (std::invalid_argument) {
rightnumber = false;
}
if (!rightnumber) {
std::cout << "Please, enter parameter " << num << ": ";
std::cin >> res;
}
}
else {
std::cout << "Please, enter parameter " << num << ": ";
std::cin >> res;
}
return res;
}
long int takeintpar (int argc, char* argv[], int num) {
long int res = 0;
if (argc > num) {
bool rightnumber = true;
try {
res = std::stoi (argv[num]);
} catch (std::invalid_argument) {
rightnumber = false;
}
if (!rightnumber) {
std::cout << "Please, enter parameter " << num << ": ";
std::cin >> res;
}
}
else {
std::cout << "Please, enter parameter " << num << ": ";
std::cin >> res;
}
return res;
}
int main (int argc, char* argv[]) {
double pars[totalpars];
pars[0] = takefloatpar(argc, argv, 1);
pars[1] = takefloatpar(argc, argv, 2);
double vars[totalvars];
vars[0] = takefloatpar(argc, argv, 3);
vars[1] = takefloatpar(argc, argv, 4);
double time = takefloatpar(argc, argv, 5);
double trans = takefloatpar(argc, argv, 6);
double timestep = takefloatpar(argc, argv, 7);
singleparset(vars,pars,2,time,trans,timestep);
return 0;
}
p. s. Чтобы сразу увидеть новый материал в моем блоге в своей ленте, подписывайтесь! Буду рад комментариям, вопросам, предложениям.
