Next: , Previous: , Up: Data processing   [Contents][Index]


7.11 Глобальные функции

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

Команда MGL: transform DAT 'type' real imag
Global function: mglData mglTransform (const mglDataA &real, const mglDataA &imag, const char *type)
Функция С: HMDT mgl_transform (HCDT real, HCDT imag, const char *type)

Выполняет интегральное преобразование комплексных данных real, imag в выбранном направлении и возвращает модуль результата. Порядок и тип преобразований задается строкой type: первый символ для x-направления, второй для y-направления, третий для z-направления. Возможные символы: ‘f’ – прямое преобразование Фурье, ‘i’ – обратное преобразование Фурье, ‘s’ – синус преобразование, ‘c’ – косинус преобразование, ‘h’ – преобразование Ханкеля, ‘n’ или ‘ ’ – нет преобразования.

Команда MGL: transforma DAT 'type' ampl phase
Global function: mglData mglTransformA const mglDataA &ampl, const mglDataA &phase, const char *type)
Функция С: HMDT mgl_transform_a HCDT ampl, HCDT phase, const char *type)

Аналогично предыдущему с заданными амплитудой ampl и фазой phase комплексных чисел.

Команда MGL: fourier reDat imDat 'dir'
Команда MGL: fourier complexDat 'dir'
Global function: void mglFourier const mglDataA &re, const mglDataA &im, const char *dir)
Функция С: void mgl_data_fourier HCDT re, HCDT im, const char *dir)

Выполняет Фурье преобразование для комплексных данных re+i*im в направлениях dir. Результат помещается обратно в массивы re и im.

Команда MGL: stfad RES real imag dn ['dir'='x']
Global function: mglData mglSTFA (const mglDataA &real, const mglDataA &imag, int dn, char dir='x')
Функция С: HMDT mgl_data_stfa (HCDT real, HCDT imag, int dn, char dir)

Выполняет оконное преобразование Фурье длиной dn для комплексных данных real, imag и возвращает модуль результата. Например, для dir=‘x’ результат будет иметь размер {int(nx/dn), dn, ny} и будет равен res[i,j,k]=|\sum_d^dn exp(I*j*d)*(real[i*dn+d,k]+I*imag[i*dn+d,k])|/dn.

Команда MGL: tridmat RES ADAT BDAT CDAT DDAT 'how'
Global function: mglData mglTridMat (const mglDataA &A, const mglDataA &B, const mglDataA &C, const mglDataA &D, const char *how)
Global function: mglDataC mglTridMatC (const mglDataA &A, const mglDataA &B, const mglDataA &C, const mglDataA &D, const char *how)
Функция С: HMDT mgl_data_tridmat (HCDT A, HCDT B, HCDT C, HCDT D, const char*how)
Функция С: HADT mgl_datac_tridmat (HCDT A, HCDT B, HCDT C, HCDT D, const char*how)

Возвращает решение трехдиагональной системы уравнений A[i]*x[i-1]+B[i]*x[i]+C[i]*x[i+1]=D[i]. Строка how может содержать:

Размеры массивов A, B, C должны быть одинаковы. Также их размерности должны совпадать со всеми или с "младшими" размерностями массива D. См. раздел PDE solving hints, для примеров кода и графика.

Команда MGL: pde RES 'ham' ini_re ini_im [dz=0.1 k0=100]
Global function: mglData mglPDE (HMGL gr, const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, mreal dz=0.1, mreal k0=100, const char *opt="")
Функция С: HMDT mgl_pde_solve (HMGL gr, const char *ham, HCDT ini_re, HCDT ini_im, mreal dz, mreal k0, const char *opt)

Решает уравнение в частных производных du/dz = i*k0*ham(p,q,x,y,z,|u|)[u], где p=-i/k0*d/dx, q=-i/k0*d/dy – псевдо-дифференциальные операторы. Параметры ini_re, ini_im задают начальное распределение поля. Координаты в уравнении и в решении полагаются в диапазоне осей координат. Замечу, что внутри этот диапазон увеличивается в 3/2 раза для уменьшения отражения от границ расчетного интервала. Параметр dz задает шаг по эволюционной координате z. В данный момент использован упрощенный алгоритм, когда все “смешанные” члена (типа ‘x*p’->x*d/dx) исключаются. Например, в 2D случае это функции типа ham = f(p,z) + g(x,z,u). При этом допускаются коммутирующие комбинации (типа ‘x*q’->x*d/dy). Переменная ‘u’ используется для обозначения амплитуды поля |u|. Это позволяет решать нелинейные задачи – например, нелинейное уравнение Шредингера ham='p^2+q^2-u^2'. Также можно указать мнимую часть для поглощения (типа ham = 'p^2+i*x*(x>0)'), но только если зависимость от ‘i’ линейная, т.е. ham = hre+i*him. См. раздел PDE solving hints, для примеров кода и графика.

Команда MGL: ray RES 'ham' x0 y0 z0 p0 q0 v0 [dt=0.1 tmax=10]
Global function: mglData mglRay (const char *ham, mglPoint r0, mglPoint p0, mreal dt=0.1, mreal tmax=10)
Функция С: HMDT mgl_ray_trace (const char *ham, mreal x0, mreal y0, mreal z0, mreal px, mreal py, mreal pz, mreal dt, mreal tmax)

Решает систему геометрооптических уравнений dr/dt = d ham/dp, dp/dt = -d ham/dr. Это гамильтоновы уравнения для траектории частицы в 3D случае. Гамильтониан ham может зависеть от координат ‘x’, ‘y’, ‘z’, импульсов ‘p’=px, ‘q’=py, ‘v’=pz и времени ‘t’: ham = H(x,y,z,p,q,v,t). Начальная точка (при t=0) задается переменными {x0, y0, z0, p0, q0, v0}. Параметры dt и tmax задают шаг и максимальное время интегрирования. Результат – массив {x,y,z,p,q,v,t} с размером {7 * int(tmax/dt+1) }.

Команда MGL: ode RES 'df' 'var' ini [dt=0.1 tmax=10]
Global function: mglData mglODE (const char *df, const char *var, const mglDataA &ini, mreal dt=0.1, mreal tmax=10)
Функция С: HMDT mgl_ode_solve_str (const char *df, const char *var, HCDT ini, mreal dt, mreal tmax)
Функция С: HMDT mgl_ode_solve (void (*df)(const mreal *x, mreal *dx, void *par), int n, const mreal *ini, mreal dt, mreal tmax)
Функция С: HMDT mgl_ode_solve_ex (void (*df)(const mreal *x, mreal *dx, void *par), int n, const mreal *ini, mreal dt, mreal tmax, void (*bord)(mreal *x, const mreal *xprev, void *par))

Решает систему обыкновенных дифференциальных уравнений dx/dt = df(x). Функции df могут быть заданны строкой с разделенными ’;’ формулами (аргумент var задает символы для переменных x[i]) или указателем на функцию, которая заполняет dx по заданным значениям x. Параметры ini, dt, tmax задают начальные значения, шаг и максимальное время интегрирования. Результат – массив с размером {n * int(tmax/dt+1)}.

Команда MGL: qo2d RES 'ham' ini_re ini_im ray [r=1 k0=100 xx yy]
Global function: mglData mglQO2d (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mreal r=1, mreal k0=100, mglData *xx=0, mglData *yy=0)
Global function: mglData mglQO2d (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mglData &xx, mglData &yy, mreal r=1, mreal k0=100)
Global function: mglDataC mglQO2dc (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mreal r=1, mreal k0=100)
Global function: mglDataC mglQO2dc (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mglData &xx, mglData &yy, mreal r=1, mreal k0=100)
Функция С: HMDT mgl_qo2d_solve (const char *ham, HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy)
Функция С: HADT mgl_qo2d_solve_c (const char *ham, HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy)
Функция С: HMDT mgl_qo2d_func (dual (*ham)(mreal u, mreal x, mreal y, mreal px, mreal py, void *par), HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy)
Функция С: HADT mgl_qo2d_func_c (dual (*ham)(mreal u, mreal x, mreal y, mreal px, mreal py, void *par), HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy)

Решает уравнение в частных производных du/dt = i*k0*ham(p,q,x,y,|u|)[u] в сопровождающей системе координат, где p=-i/k0*d/dx, q=-i/k0*d/dy – псевдо-дифференциальные операторы. Параметры ini_re, ini_im задают начальное распределение поля. Параметр ray задает опорный луч для сопровождающей системы координат. Можно использовать луч найденный с помощью ray. Опорный луч должен быть достаточно гладкий, чтобы система координат была однозначной и для исключения ошибок интегрирования. Если массивы xx и yy указаны, то в них записываются декартовы координаты для каждой точки найденного решения. См. также pde, qo3d. См. раздел PDE solving hints, для примеров кода и графика.

Команда MGL: qo3d RES 'ham' ini_re ini_im ray [r=1 k0=100 xx yy zz]
Global function: mglData mglQO3d (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mreal r=1, mreal k0=100)
Global function: mglData mglQO3d (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mglData &xx, mglData &yy, mglData &zz, mreal r=1, mreal k0=100)
Global function: mglDataC mglQO3dc (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mreal r=1, mreal k0=100)
Global function: mglDataC mglQO3dc (const char *ham, const mglDataA &ini_re, const mglDataA &ini_im, const mglDataA &ray, mglData &xx, mglData &yy, mglData &zz, mreal r=1, mreal k0=100)
Функция С: HMDT mgl_qo3d_solve (const char *ham, HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy, HMDT zz)
Функция С: HADT mgl_qo3d_solve_c (const char *ham, HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy, HMDT zz)
Функция С: HMDT mgl_qo3d_func (dual (*ham)(mreal u, mreal x, mreal y, mreal z, mreal px, mreal py, mreal pz, void *par), HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy, HMDT zz)
Функция С: HADT mgl_qo3d_func_c (dual (*ham)(mreal u, mreal x, mreal y, mreal z, mreal px, mreal py, mreal pz, void *par), HCDT ini_re, HCDT ini_im, HCDT ray, mreal r, mreal k0, HMDT xx, HMDT yy, HMDT zz)

Решает уравнение в частных производных du/dt = i*k0*ham(p,q,v,x,y,z,|u|)[u] в сопровождающей системе координат, где p=-i/k0*d/dx, q=-i/k0*d/dy, v=-i/k0*d/dz – псевдо-дифференциальные операторы. Параметры ini_re, ini_im задают начальное распределение поля. Параметр ray задает опорный луч для сопровождающей системы координат. Можно использовать луч найденный с помощью ray. Опорный луч должен быть достаточно гладкий, чтобы система координат была однозначной и для исключения ошибок интегрирования. Если массивы xx, yy и zz указаны, то в них записываются декартовы координаты для каждой точки найденного решения. См. также pde, qo2d.

Команда MGL: jacobian RES xdat ydat [zdat]
Global function: mglData mglJacobian (const mglDataA &x, const mglDataA &y)
Global function: mglData mglJacobian (const mglDataA &x, const mglDataA &y, const mglDataA &z)
Функция С: HMDT mgl_jacobian_2d (HCDT x, HCDT y)
Функция С: HMDT mgl_jacobian_3d (HCDT x, HCDT y, HCDT z)

Вычисляет якобиан преобразования {i,j,k} в {x,y,z}, где координаты {i,j,k} полагаются нормированными в интервал [0,1]. Якобиан находится по формуле det||dr_\alpha/d\xi_\beta||, где r={x,y,z} и \xi={i,j,k}. Все размерности всех массивов должны быть одинаковы. Данные должны быть трехмерными если указаны все 3 массива {x,y,z} или двумерными если только 2 массива {x,y}.

Команда MGL: triangulation RES xdat ydat [zdat]
Global function: mglData mglTriangulation (const mglDataA &x, const mglDataA &y)
Global function: mglData mglTriangulation (const mglDataA &x, const mglDataA &y, const mglDataA &z)
Функция С: HMDT mgl_triangulation_2d (HCDT x, HCDT y)
Функция С: HMDT mgl_triangulation_3d (HCDT x, HCDT y, HCDT z)

Выполняет триангуляцию для произвольно расположенных точек с координатами {x,y,z} (т.е. находит треугольники, соединяющие точки). Первая размерность всех массивов должна быть одинакова x.nx=y.nx=z.nx. Получившийся массив можно использовать в triplot или tricont для визуализации реконструированной поверхности. См. раздел Making regular data, для примеров кода и графика.

Global function: mglData mglGSplineInit (const mglDataA &x, const mglDataA &y)
Global function: mglDataC mglGSplineCInit (const mglDataA &x, const mglDataA &y)
Функция С: HMDT mgl_gspline_init (HCDT x, HCDT y)
Функция С: HADT mgl_gsplinec_init (HCDT x, HCDT y)

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

Global function: mreal mglGSpline (const mglDataA &coef, mreal dx, mreal *d1=0, mreal *d2=0)
Global function: dual mglGSplineC (const mglDataA &coef, mreal dx, dual *d1=0, dual *d2=0)
Функция С: mreal mgl_gspline (HCDT coef, mreal dx, mreal *d1, mreal *d2)
Функция С: dual mgl_gsplinec (HCDT coef, mreal dx, dual *d1, dual *d2)

Вычисляет глобальный кубический сплайн (а также 1ую и 2ую производные d1, d2 если они не NULL), используя коэффициенты coef в точке dx+x0 (здесь x0 – 1ый элемент массива x в функции mglGSpline*Init()).

Команда MGL: ifs2d RES dat num [skip=20]
Global function: mglData mglIFS2d (const mglDataA &dat, long num, long skip=20)
Функция С: HMDT mgl_data_ifs_2d (HCDT dat, long num, long skip)

Находит num точек {x[i]=res[0,i], y[i]=res[1,i]} фрактала с использованием итерационной системы функций (IFS). Матрица dat используется для генерации в соответствии с формулами

x[i+1] = dat[0,i]*x[i] + dat[1,i]*y[i] + dat[4,i];
y[i+1] = dat[2,i]*x[i] + dat[3,i]*y[i] + dat[5,i];

Значение dat[6,i] – весовой коэффициент для i-ой строки матрицы dat. Первые skip итераций будут опущены. Массив dat должен иметь размер по x больше или равный 7. См. раздел IFS sample, для примеров кода и графика.

Команда MGL: ifs3d RES dat num [skip=20]
Global function: mglData mglIFS3d (const mglDataA &dat, long num, long skip=20)
Функция С: HMDT mgl_data_ifs_3d (HCDT dat, long num, long skip)

Находит num точек {x[i]=res[0,i], y[i]=res[1,i], z[i]=res[2,i]} фрактала с использованием итерационной системы функций (IFS). Матрица dat используется для генерации в соответствии с формулами

x[i+1] = dat[0,i]*x[i] + dat[1,i]*y[i] + dat[2,i]*z[i] + dat[9,i];
y[i+1] = dat[3,i]*x[i] + dat[4,i]*y[i] + dat[5,i]*z[i] + dat[10,i];
z[i+1] = dat[6,i]*x[i] + dat[7,i]*y[i] + dat[8,i]*z[i] + dat[10,i];

Значение dat[6,i] – весовой коэффициент для i-ой строки матрицы dat. Первые skip итераций будут опущены. Массив dat должен иметь размер по x больше или равный 13. См. раздел IFS sample, для примеров кода и графика.

Команда MGL: ifsfile RES 'fname' 'name' num [skip=20]
Global function: mglData mglIFSfile (const char *fname, const char *name, long num, long skip=20)
Функция С: HMDT mgl_data_ifs_file (const char *fname, const char *name, long num, long skip)

Считывает параметры фрактала name из файла fname и находит num точек для него. Первые skip итераций будут опущены. См. также ifs2d, ifs3d.

Файл IFS может содержать несколько записей. Каждая запись содержит имя фрактала (‘binary’ в примере ниже) и тело в фигурных скобках {} с параметрами фрактала. Символ ‘;’ начинает комментарий. Если имя содержит ‘(3D)’ или ‘(3d)’, то определен 3d IFS фрактал. Пример содержит два фрактала: ‘binary’ – обычный 2d фрактал, и ‘3dfern (3D)’ – 3d фрактал.

 binary
 { ; comment allowed here
  ; and here
  .5  .0 .0 .5 -2.563477 -0.000003 .333333   ; also comment allowed here
  .5  .0 .0 .5  2.436544 -0.000003 .333333
  .0 -.5 .5 .0  4.873085  7.563492 .333333
  }

 3dfern (3D) {
   .00  .00 0 .0 .18 .0 0  0.0 0.00 0 0.0 0 .01
   .85  .00 0 .0 .85 .1 0 -0.1 0.85 0 1.6 0 .85
   .20 -.20 0 .2 .20 .0 0  0.0 0.30 0 0.8 0 .07
  -.20  .20 0 .2 .20 .0 0  0.0 0.30 0 0.8 0 .07
  }




Next: , Previous: , Up: Data processing   [Contents][Index]