Эти функции не методы класса mglData
, но они дают дополнительные возможности по обработке данных. Поэтому я поместил их в эту главу.
DAT 'type' real imag
¶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’ или ‘ ’ – нет преобразования.
DAT 'type' ampl phase
¶mglData
mglTransformA const mglDataA &
ampl, const mglDataA &
phase, const char *
type)
¶HMDT
mgl_transform_a HCDT
ampl, HCDT
phase, const char *
type)
¶Аналогично предыдущему с заданными амплитудой ampl и фазой phase комплексных чисел.
reDat imDat 'dir'
¶complexDat 'dir'
¶void
mglFourier const mglDataA &
re, const mglDataA &
im, const char *
dir)
¶mglDataC
: void
FFT (const char *
dir)
¶void
mgl_data_fourier HCDT
re, HCDT
im, const char *
dir)
¶void
mgl_datac_fft (HADT
dat, const char *
dir)
¶Выполняет Фурье преобразование для комплексных данных re+i*im в направлениях dir. Результат помещается обратно в массивы re и im. Если dir содержит ‘i’, то выполняется обратное преобразование Фурье.
RES real imag dn
['dir'='x']
¶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.
dat xdat ydat
¶mglData
mglTriangulation (const mglDataA &
x, const mglDataA &
y)
¶void
mgl_triangulation_2d (HCDT
x, HCDT
y)
¶Выполняет триангуляцию Делоне для точек на плоскости и возвращает массив, пригодный для triplot и tricont. См. раздел Making regular data, для примеров кода и графика.
RES ADAT BDAT CDAT DDAT 'how'
¶mglData
mglTridMat (const mglDataA &
A, const mglDataA &
B, const mglDataA &
C, const mglDataA &
D, const char *
how)
¶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, для примеров кода и графика.
RES 'ham' ini_re ini_im [dz=0.1 k0=100
]
¶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=""
)
¶mglDataC
mglPDEc (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)
¶HADT
mgl_pde_solve_c (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)'
). См. также apde, qo2d, qo3d. См. раздел PDE solving hints, для примеров кода и графика.
RES 'ham' ini_re ini_im [dz=0.1 k0=100
]
¶mglData
mglAPDE (HMGL
gr, const char *
ham, const mglDataA &
ini_re, const mglDataA &
ini_im, mreal
dz=0.1
, mreal
k0=100
, const char *
opt=""
)
¶mglDataC
mglAPDEc (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_adv (HMGL
gr, const char *
ham, HCDT
ini_re, HCDT
ini_im, mreal
dz, mreal
k0, const char *
opt)
¶HADT
mgl_pde_solve_adv_c (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. Используется достаточно сложный и медленный алгоритм, способный учесть одновременное влияние пространственной дисперсии и неоднородности среды [см. А.А. Балакин, Е.Д. Господчиков, А.Г. Шалашов, Письма ЖЭТФ 104, 701 (2016)]. Переменная ‘u’ используется для обозначения амплитуды поля |u|. Это позволяет решать нелинейные задачи – например, нелинейное уравнение Шредингера ham='p^2+q^2-u^2'
. Также можно указать мнимую часть для поглощения (типа ham = 'p^2+i*x*(x>0)'
). См. также apde. См. раздел PDE solving hints, для примеров кода и графика.
RES 'ham' x0 y0 z0 p0 q0 v0 [dt=0.1 tmax=10]
¶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) }.
RES 'df' 'var' ini [dt=0.1 tmax=10
'jump'='']
¶mglData
mglODE (const char *
df, const char *
var, const mglDataA &
ini, mreal
dt=0.1
, mreal
tmax=10
, const char *
jump=""
)
¶mglDataC
mglODEc (const char *
df, const char *
var, const mglDataA &
ini, mreal
dt=0.1
, mreal
tmax=10
, const char *
jump=""
)
¶HMDT
mgl_ode_solve_str (const char *
df, const char *
var, HCDT
ini, mreal
dt, mreal
tmax)
¶HMDT
mgl_ode_solve_str_b (const char *
df, const char *
var, HCDT
ini, mreal
dt, mreal
tmax, const char *
jump)
¶HADT
mgl_ode_solve_str_c (const char *
df, const char *
var, HCDT
ini, mreal
dt, mreal
tmax)
¶HADT
mgl_ode_solve_str_cb (const char *
df, const char *
var, HCDT
ini, mreal
dt, mreal
tmax, const char *
jump)
¶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 (*
jump)(mreal *x, const mreal *xprev, void *par)
)
¶Решает систему обыкновенных дифференциальных уравнений dx/dt = df(x). Функции df могут быть заданны строкой с разделенными ’;’ формулами (аргумент var задает символы для переменных x[i]) или указателем на функцию, которая заполняет dx
по заданным значениям x
. Параметры ini, dt, tmax задают начальные значения, шаг и максимальное время интегрирования. Функция обрывает расчет при появлении значений NAN
или INF
. Результат – массив размером {n * Nt}, где Nt <= int(tmax/dt+1).
Если dt*tmax<0, то включается регуляризация, при которой уравнение заменяется на dx/ds = df(x)/max(|df(x)|) для более аккуратного прохождения областей резкого изменения df и более быстрого расчета в области малых df. Здесь s – новое "время". При этом реальное время определяется уравнением dt/ds=max(|df(x)|). Поэтому для вывода реального времени его следует явно добавить в уравнение, например ‘ode res 'y;-sin(x);1' 'xyt' [3,0] 0.3 -100’. Оно также сохранит точность вблизи стационарных точек (т.е. при малых df в периодическом движении).
Функция jump позволяет задать особое поведение, например, отражение от границы или прыжки из-за толчков.
RES 'df' 'var' 'brd' ini [dt=0.1 tmax=10
]
¶mglData
mglODEs (const char *
df, const char *
var, char
brd, const mglDataA &
ini, mreal
dt=0.1
, mreal
tmax=10
)
¶mglDataC
mglODEcs (const char *
df, const char *
var, char
brd, const mglDataA &
ini, mreal
dt=0.1
, mreal
tmax=10
)
¶HMDT
mgl_ode_solve_set (const char *
df, const char *
var, char
brd, HCDT
ini, mreal
dt, mreal
tmax)
¶HADT
mgl_ode_solve_set_c (const char *
df, const char *
var, char
brd, HCDT
ini, mreal
dt, mreal
tmax)
¶Решает разностную аппроксимацию уравнений в частных производных как систему обыкновенных дифференциальных уравнений dx/dt = df(x,j). Функции df могут быть заданны строкой с разделенными ’;’ формулами, зависящими от индекса j и текущего времени ‘t’. Аргумент var задает символы для переменных x[i]. Параметр brd задает тип граничных условий по j: ‘0’ или ‘z’ – нулевые, ‘1’ или ‘c’ – постоянные, ‘2’ или ‘l’ – линейные (ноль лапласиана), ‘3’ или ‘s’ – квадратичные, ‘4’ или ‘e’ – экспоненциальные, ‘5’ или ‘g’ – гауссовы. Последние два типа условий применимы только в случае комплексных переменных. Параметры ini, dt, tmax задают начальные значения, шаг и максимальное время интегрирования. Функция обрывает расчет при появлении значений NAN
или INF
. Результат – массив размером {n * Nt}, где Nt <= int(tmax/dt+1). Например, разностная аппроксимация уравнения диффузии с нулевыми граничными условиями находится вызовом: ‘ode res 'u(j+1)-2*u(j)+u(j-1)' 'u' '0' u0’, где ‘u0’ – массив начальных значений.
Если dt*tmax<0, то включается регуляризация, при которой уравнение заменяется на dx/ds = df(x)/max(|df(x)|) для более аккуратного прохождения областей резкого изменения df и более быстрого расчета в области малых df. Здесь s – новое "время". При этом реальное время определяется уравнением dt/ds=max(|df(x)|). Поэтому для вывода реального времени его следует явно добавить в уравнение, например ‘ode res 'y;-sin(x);1' 'xyt' [3,0] 0.3 -100’. Оно также сохранит точность вблизи стационарных точек (т.е. при малых df в периодическом движении).
RES 'ham' ini_re ini_im ray [r=1 k0=100
xx yy]
¶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
)
¶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
)
¶mglDataC
mglQO2dc (const char *
ham, const mglDataA &
ini_re, const mglDataA &
ini_im, const mglDataA &
ray, mreal
r=1
, mreal
k0=100
)
¶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, для примеров кода и графика.
RES 'ham' ini_re ini_im ray [r=1 k0=100
xx yy zz]
¶mglData
mglQO3d (const char *
ham, const mglDataA &
ini_re, const mglDataA &
ini_im, const mglDataA &
ray, mreal
r=1
, mreal
k0=100
)
¶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
)
¶mglDataC
mglQO3dc (const char *
ham, const mglDataA &
ini_re, const mglDataA &
ini_im, const mglDataA &
ray, mreal
r=1
, mreal
k0=100
)
¶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.
RES xdat ydat [zdat]
¶mglData
mglJacobian (const mglDataA &
x, const mglDataA &
y)
¶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}.
RES xdat ydat [zdat]
¶mglData
mglTriangulation (const mglDataA &
x, const mglDataA &
y)
¶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, для примеров кода и графика.
mglData
mglGSplineInit (const mglDataA &
x, const mglDataA &
y)
¶mglDataC
mglGSplineCInit (const mglDataA &
x, const mglDataA &
y)
¶HMDT
mgl_gspline_init (HCDT
x, HCDT
y)
¶HADT
mgl_gsplinec_init (HCDT
x, HCDT
y)
¶Подготавливает коэффициенты для глобального кубического сплайна.
mreal
mglGSpline (const mglDataA &
coef, mreal
dx, mreal *
d1=0
, mreal *
d2=0
)
¶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()
).
RES dat num
[skip=20
]
¶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. См. также ifs3d, flame2d. См. раздел Sample ‘ifs2d’, для примеров кода и графика.
RES dat num
[skip=20
]
¶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[11,i];
Значение dat[12,i]
– весовой коэффициент для i-ой строки матрицы dat. Первые skip итераций будут опущены. Массив dat должен иметь размер по x больше или равный 13. См. также ifs2d. См. раздел Sample ‘ifs3d’, для примеров кода и графика.
RES 'fname' 'name' num
[skip=20
]
¶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 фрактал. См. также ifs2d, ifs3d.
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 }
RES dat func num
[skip=20
]
¶mglData
mglFlame2d (const mglDataA &
dat, const mglDataA &
func, long
num, long
skip=20
)
¶HMDT
mgl_data_flame_2d (HCDT
dat, HCDT
func, long
num, long
skip)
¶Находит num точек {x[i]=res[0,i], y[i]=res[1,i]} фрактала с использованием итерационной системы функций (IFS). Массив func задает идентификатор функции (func[0,i,j]), ее вес (func[0,i,j]) и аргументы (func[2 ... 5,i,j]). Матрица dat используется для преобразования координат для аргументов функции. Результирующее преобразование имеет вид:
xx = dat[0,i]*x[j] + dat[1,j]*y[i] + dat[4,j]; yy = dat[2,i]*x[j] + dat[3,j]*y[i] + dat[5,j]; x[j+1] = sum_i @var{func}[1,i,j]*@var{func}[0,i,j]_x(xx, yy; @var{func}[2,i,j],...,@var{func}[5,i,j]); y[j+1] = sum_i @var{func}[1,i,j]*@var{func}[0,i,j]_y(xx, yy; @var{func}[2,i,j],...,@var{func}[5,i,j]);
Значение dat[6,i]
– весовой коэффициент для i-ой строки матрицы dat. Первые skip итераций будут опущены. Массив dat должен иметь размер по x больше или равный 7.
Доступные идентификаторы функций: mglFlame2d_linear=0, mglFlame2d_sinusoidal, mglFlame2d_spherical, mglFlame2d_swirl, mglFlame2d_horseshoe,
mglFlame2d_polar, mglFlame2d_handkerchief,mglFlame2d_heart, mglFlame2d_disc, mglFlame2d_spiral,
mglFlame2d_hyperbolic, mglFlame2d_diamond, mglFlame2d_ex, mglFlame2d_julia, mglFlame2d_bent,
mglFlame2d_waves, mglFlame2d_fisheye, mglFlame2d_popcorn, mglFlame2d_exponential, mglFlame2d_power,
mglFlame2d_cosine, mglFlame2d_rings, mglFlame2d_fan, mglFlame2d_blob, mglFlame2d_pdj,
mglFlame2d_fan2, mglFlame2d_rings2, mglFlame2d_eyefish, mglFlame2d_bubble, mglFlame2d_cylinder,
mglFlame2d_perspective, mglFlame2d_noise, mglFlame2d_juliaN, mglFlame2d_juliaScope, mglFlame2d_blur,
mglFlame2d_gaussian, mglFlame2d_radialBlur, mglFlame2d_pie, mglFlame2d_ngon, mglFlame2d_curl,
mglFlame2d_rectangles, mglFlame2d_arch, mglFlame2d_tangent, mglFlame2d_square, mglFlame2d_blade,
mglFlame2d_secant, mglFlame2d_rays, mglFlame2d_twintrian, mglFlame2d_cross, mglFlame2d_disc2,
mglFlame2d_supershape, mglFlame2d_flower, mglFlame2d_conic, mglFlame2d_parabola, mglFlame2d_bent2,
mglFlame2d_bipolar, mglFlame2d_boarders, mglFlame2d_butterfly, mglFlame2d_cell, mglFlame2d_cpow,
mglFlame2d_curve, mglFlame2d_edisc, mglFlame2d_elliptic, mglFlame2d_escher, mglFlame2d_foci,
mglFlame2d_lazySusan, mglFlame2d_loonie, mglFlame2d_preBlur, mglFlame2d_modulus, mglFlame2d_oscope,
mglFlame2d_polar2, mglFlame2d_popcorn2, mglFlame2d_scry, mglFlame2d_separation, mglFlame2d_split,
mglFlame2d_splits, mglFlame2d_stripes, mglFlame2d_wedge, mglFlame2d_wedgeJulia, mglFlame2d_wedgeSph,
mglFlame2d_whorl, mglFlame2d_waves2, mglFlame2d_exp, mglFlame2d_log, mglFlame2d_sin,
mglFlame2d_cos, mglFlame2d_tan, mglFlame2d_sec, mglFlame2d_csc, mglFlame2d_cot,
mglFlame2d_sinh, mglFlame2d_cosh, mglFlame2d_tanh, mglFlame2d_sech, mglFlame2d_csch,
mglFlame2d_coth, mglFlame2d_auger, mglFlame2d_flux.
Значение dat[6,i]
– весовой коэффициент для i-ой строки матрицы dat. Первые skip итераций будут опущены. Размеры массивов должны удовлетворять требованиям: dat.nx>=7, func.nx>=2 и func.nz=dat.ny. См. также ifs2d, ifs3d. См. раздел Sample ‘flame2d’, для примеров кода и графика.