NURBS: Где я могу найти эти две функции полезности линейной алгебры? - PullRequest
0 голосов
/ 16 мая 2018

Я работаю через Книгу NURBS Пигла и Тиллера. Для алгоритма глобальной интерполяции требуется, чтобы вы предоставили две служебные программы для решения системы линейных уравнений:

LUDecomposition(A, q, sbw) для разложения матрицы коэффициентов q x q с шириной полосы пропускания sbw на нижнюю и верхнюю треугольные составляющие; для простоты мы предполагаем, что A является квадратным массивом q x q, но используется утилита, которая хранит только ненулевую полосу.

ForwardBackward(A, q, sbw, rhs, sol) для выполнения прямого / обратного замещения (см. [Press88]); rhs [] - правая часть системы (координаты Q_k), а sol [] - вектор решения (координаты P_i).

Проверяя ссылку Press88, я обнаружил, что это Числовые рецепты на C . Я должен быть в состоянии переделать алгоритм в этой книге, чтобы получить функцию ForwardBackward, но что касается LUDecomposition, где я могу найти тот, который работает для особого случая матрицы с диагональными полосами?

1 Ответ

0 голосов
/ 02 июня 2018

Код Matlab для LU decomp Трехугольная матрица

 function [u1,d1,l1] = decomt(u,d,l)
    %length of diagonal
    n=length(d);
    u1=u;
    d1 = d;
    l1=l;

    %perform LU decomp

    d1(1) = d(1);
    for i =2:n
        l1(i-1) = l(i-1)/d1(i-1); %update lower triangle
        d1(i)= d(i) - (l(i-1)/d1(i-1))*u(i-1); % update diagonal
    end
    end

Чтобы решить, передний саб и задний саб.

function [x] = solvet(u,d,l,b) 

n=length(d);
x = (1:n);
y =(1:n);

% Solve tridiag LUx=b

% Step 1 Solve Ly=b for y

y(1) = b(1)

for i=2:n
    y(i) = b(i) - l(i-1)*y(i-1);
end

%Step 2 : Solve Ux=y for x
x(n) = y(n)/d(n);
for i=(n-1):-1:1
    x(i) = (y(i)-u(i)*x(i+1))/d(i);
end
end
...