changeset 660:184833fe4c0e feature/grids

Add monomial, vandermonde and stencil cretation tools
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 14 Nov 2017 16:46:07 +0100
parents 11a39b274260
children 905612b7f3d4
files mononomial.m stencilEquation.m vandermonde.m
diffstat 3 files changed, 31 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mononomial.m	Tue Nov 14 16:46:07 2017 +0100
@@ -0,0 +1,8 @@
+function y = mononomial(x, k)
+    if k < 0
+        y = x*0;
+        return
+    end
+    y = x.^k/factorial(k);
+end
+
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/stencilEquation.m	Tue Nov 14 16:46:07 2017 +0100
@@ -0,0 +1,13 @@
+% Find the equation for the stencil for d^k/dx^k
+function [A,b] = stencilEquation(k, offsets, order)
+    q = sym('q', [1, length(offsets)]);
+
+    p = 0:(order-1+k);
+
+    v     = vandermonde(offsets, p);
+    vdiff = vandermonde(      0, p-k);
+
+    eq = q*v == vdiff;
+
+    [A,b] = equationsToMatrix(eq, q);
+end
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/vandermonde.m	Tue Nov 14 16:46:07 2017 +0100
@@ -0,0 +1,10 @@
+% Create vandermonde matrix for points x and polynomials of order p
+% x and p are vectors
+% v is a length(x) by length(p) matrix
+function V = vandermonde(x, p)
+    V = sym(zeros(length(x), length(p))); % Is there a way to make this work for both double and sym
+
+    for i = 1:length(p)
+        V(:, i) = mononomial(x,p(i));
+    end
+end