34 files changed, 384 insertions, 6 deletions
diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..9d0b71a
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,2 @@
+build
+dist
diff --git a/Makefile b/Makefile
new file mode 100644
index 0000000..a76a59f
--- /dev/null
+++ b/Makefile
@@ -0,0 +1,41 @@
+TEST_SRC := test/main.c
+SRC_DIR := src
+OBJ_DIR := build
+BIN_DIR := dist
+LIB_DIR := lib
+
+TEST_EXE := $(BIN_DIR)/lizfcm.test
+EXE := $(BIN_DIR)/lizfcm
+LIBRARY := $(LIB_DIR)/lizfcm.a
+SRC := $(wildcard $(SRC_DIR)/*.c)
+OBJ := $(SRC:$(SRC_DIR)/%.c=$(OBJ_DIR)/%.o)
+
+CPPFLAGS := -Iinc -MMD -MP
+CFLAGS   := -Wall
+LDFLAGS  := 
+LDLIBS   := -lm
+
+.PHONY: all clean
+
+all: $(TEST_EXE)
+
+$(TEST_EXE): $(LIBRARY)
+	$(CC) $(CPPFLAGS) $(CFLAGS) $(LIBRARY) $(TEST_SRC) -o $@
+
+$(LIBRARY): $(EXE)
+	ar rcs $(LIBRARY) $(OBJ_DIR)/*.o
+	ranlib $(LIBRARY)
+
+$(EXE): $(OBJ) | $(BIN_DIR)
+	$(CC) $(LDFLAGS) $^ $(LDLIBS) -o $@
+
+$(OBJ_DIR)/%.o: $(SRC_DIR)/%.c | $(OBJ_DIR)
+	$(CC) $(CPPFLAGS) $(CFLAGS) -c $< -o $@
+
+$(BIN_DIR) $(OBJ_DIR) $(LIB_DIR):
+	mkdir -p $@
+
+clean:
+	@$(RM) -r $(BIN_DIR) $(OBJ_DIR)
+
+-include $(OBJ:.o=.d)
diff --git a/compile_archive.sh b/compile_archive.sh
deleted file mode 100755
index f174fbe..0000000
--- a/compile_archive.sh
+++ /dev/null
@@ -1,3 +0,0 @@
-#!/bin/sh
-
-ar rcs lizfcm.a src/*
diff --git a/deprecated-cl/.main.lisp.swp b/deprecated-cl/.main.lisp.swp
new file mode 100644
index 0000000..e0a098e
--- /dev/null
+++ b/deprecated-cl/.main.lisp.swp
diff --git a/src/approx,derivative.lisp b/deprecated-cl/approx,derivative.lisp
index 631a5c0..631a5c0 100644
--- a/src/approx,derivative.lisp
+++ b/deprecated-cl/approx,derivative.lisp
diff --git a/src/approx,maceps.lisp b/deprecated-cl/approx,maceps.lisp
index e2738e4..e2738e4 100644
--- a/src/approx,maceps.lisp
+++ b/deprecated-cl/approx,maceps.lisp
diff --git a/src/approx,package.lisp b/deprecated-cl/approx,package.lisp
index a0eac80..a0eac80 100644
--- a/src/approx,package.lisp
+++ b/deprecated-cl/approx,package.lisp
diff --git a/src/lizfcm.asd b/deprecated-cl/lizfcm.asd
index 0096257..dea3ddd 100644
--- a/src/lizfcm.asd
+++ b/deprecated-cl/lizfcm.asd
@@ -9,7 +9,7 @@
                  (:file "approx,maceps" :depends-on ("approx,package"))
                  (:file "approx,derivative" :depends-on ("approx,package"))
                  (:file "approx,package")
-                 (:file "vector,least-squares")
+                 (:file "vector,least-squares" :depends-on ("vector,package"))
                  (:file "vector,distance" :depends-on ("vector,norm" "vector,package"))
                  (:file "vector,norm" :depends-on ("vector,package"))
                  (:file "vector,package")))
diff --git a/src/main.lisp b/deprecated-cl/main.lisp
index 7a8b455..7a8b455 100644
--- a/src/main.lisp
+++ b/deprecated-cl/main.lisp
diff --git a/src/tests,approx.lisp b/deprecated-cl/tests,approx.lisp
index 678ff8c..678ff8c 100644
--- a/src/tests,approx.lisp
+++ b/deprecated-cl/tests,approx.lisp
diff --git a/src/tests,maceps.lisp b/deprecated-cl/tests,maceps.lisp
index cd5ced9..cd5ced9 100644
--- a/src/tests,maceps.lisp
+++ b/deprecated-cl/tests,maceps.lisp
diff --git a/src/tests,suite.lisp b/deprecated-cl/tests,suite.lisp
index e23cfaf..e23cfaf 100644
--- a/src/tests,suite.lisp
+++ b/deprecated-cl/tests,suite.lisp
diff --git a/src/tests,table.lisp b/deprecated-cl/tests,table.lisp
index 33d4e86..33d4e86 100644
--- a/src/tests,table.lisp
+++ b/deprecated-cl/tests,table.lisp
diff --git a/src/tests,vector.lisp b/deprecated-cl/tests,vector.lisp
index 6edb1ac..6edb1ac 100644
--- a/src/tests,vector.lisp
+++ b/deprecated-cl/tests,vector.lisp
diff --git a/src/utils,package.lisp b/deprecated-cl/utils,package.lisp
index bdd5589..bdd5589 100644
--- a/src/utils,package.lisp
+++ b/deprecated-cl/utils,package.lisp
diff --git a/src/utils,table.lisp b/deprecated-cl/utils,table.lisp
index e96f37b..e96f37b 100644
--- a/src/utils,table.lisp
+++ b/deprecated-cl/utils,table.lisp
diff --git a/src/utils,within-range.lisp b/deprecated-cl/utils,within-range.lisp
index 9a0b762..9a0b762 100644
--- a/src/utils,within-range.lisp
+++ b/deprecated-cl/utils,within-range.lisp
diff --git a/src/vector,distance.lisp b/deprecated-cl/vector,distance.lisp
index 74631ce..74631ce 100644
--- a/src/vector,distance.lisp
+++ b/deprecated-cl/vector,distance.lisp
diff --git a/src/vector,least-squares.lisp b/deprecated-cl/vector,least-squares.lisp
index 687af32..687af32 100644
--- a/src/vector,least-squares.lisp
+++ b/deprecated-cl/vector,least-squares.lisp
diff --git a/src/vector,norm.lisp b/deprecated-cl/vector,norm.lisp
index aa51bce..aa51bce 100644
--- a/src/vector,norm.lisp
+++ b/deprecated-cl/vector,norm.lisp
diff --git a/src/vector,package.lisp b/deprecated-cl/vector,package.lisp
index f491908..f491908 100644
--- a/src/vector,package.lisp
+++ b/deprecated-cl/vector,package.lisp
diff --git a/inc/lizfcm.h b/inc/lizfcm.h
new file mode 100644
index 0000000..b390394
--- /dev/null
+++ b/inc/lizfcm.h
@@ -0,0 +1,30 @@
+#include "macros.h"
+#include "types.h"
+
+#ifndef LIZFCM_H
+#define LIZFCM_H
+
+extern float smaceps();
+extern double dmaceps();
+
+extern double central_derivative_at(double (*f)(double), double a, double h);
+extern double forward_derivative_at(double (*f)(double), double a, double h);
+extern double backward_derivative_at(double (*f)(double), double a, double h);
+
+extern double sum_v(Array_double *v);
+extern Array_double *add_v(Array_double *v1, Array_double *v2);
+extern Array_double *minus_v(Array_double *v1, Array_double *v2);
+extern double dot_v(Array_double *v1, Array_double *v2);
+extern double l2_norm(Array_double *v);
+extern double l1_norm(Array_double *v);
+extern double linf_norm(Array_double *v);
+
+extern double l2_distance(Array_double *v1, Array_double *v2);
+extern double l1_distance(Array_double *v1, Array_double *v2);
+extern double linf_distance(Array_double *v1, Array_double *v2);
+
+extern void format_vector_into(Array_double *v, char *s);
+
+extern Line *least_squares_lin_reg(Array_double *x, Array_double *y);
+
+#endif // LIZFCM_H
diff --git a/inc/macros.h b/inc/macros.h
new file mode 100644
index 0000000..841084d
--- /dev/null
+++ b/inc/macros.h
@@ -0,0 +1,41 @@
+#include <stddef.h>
+#include <stdlib.h>
+#include <string.h>
+
+#ifndef MACROS_H
+#define MACROS_H
+
+#define DEFINE_ARRAY(TYPE)                                                     \
+  typedef struct {                                                             \
+    TYPE *data;                                                                \
+    size_t size;                                                               \
+  } Array_##TYPE
+
+#define InitArray(TYPE, ...)                                                   \
+  ({                                                                           \
+    TYPE temp[] = __VA_ARGS__;                                                 \
+    Array_##TYPE *arr = malloc(sizeof(Array_##TYPE));                          \
+    arr->size = sizeof(temp) / sizeof(temp[0]);                                \
+    arr->data = malloc(arr->size * sizeof(TYPE));                              \
+    if (arr->data) {                                                           \
+      memcpy(arr->data, temp, arr->size * sizeof(TYPE));                       \
+    }                                                                          \
+    arr;                                                                       \
+  })
+
+#define InitArrayWithSize(TYPE, SIZE, INIT_VALUE)                              \
+  ({                                                                           \
+    Array_##TYPE *arr = malloc(sizeof(Array_##TYPE));                          \
+    arr->size = SIZE;                                                          \
+    arr->data = malloc(arr->size * sizeof(TYPE));                              \
+    if (arr->data) {                                                           \
+      for (size_t i = 0; i < arr->size; i++)                                   \
+        arr->data[i] = INIT_VALUE;                                             \
+    }                                                                          \
+    arr;                                                                       \
+  })
+
+#define c_max(x, y) (((x) >= (y)) ? (x) : (y))
+#define c_min(x, y) (((x) <= (y)) ? (x) : (y))
+
+#endif // MACROS_H
diff --git a/inc/types.h b/inc/types.h
new file mode 100644
index 0000000..4af17db
--- /dev/null
+++ b/inc/types.h
@@ -0,0 +1,18 @@
+#include "macros.h"
+#include <stddef.h>
+
+#ifndef TYPES_H
+#define TYPES_H
+
+DEFINE_ARRAY(int);
+DEFINE_ARRAY(uint32_t);
+DEFINE_ARRAY(int32_t);
+DEFINE_ARRAY(float);
+DEFINE_ARRAY(double);
+
+typedef struct Line {
+  double m;
+  double a;
+} Line;
+
+#endif // TYPES_H
diff --git a/lib/lizfcm.a b/lib/lizfcm.a
new file mode 100644
index 0000000..fe9a6ca
--- /dev/null
+++ b/lib/lizfcm.a
diff --git a/lizfcm.a b/lizfcm.a
deleted file mode 100644
index cb34117..0000000
--- a/lizfcm.a
+++ /dev/null
diff --git a/notes/Oct-11.org b/notes/Oct-11.org
new file mode 100644
index 0000000..575ea74
--- /dev/null
+++ b/notes/Oct-11.org
@@ -0,0 +1,15 @@
+* Diagonal Dominance
+Suppose that A \in R^{n \times n} is diagonally dominant then Gaussian eliminiation of A produces no zero pivot
+elements.
+
+Def. A \in R^{n \times n} is diagonally dominant if for each i=1,2,...n |a_{i,i}| \geq \Sigma_{j=1}^n |a_i,j|
+
+
+* To test solution code:
+         [[1]
+          [1]
+Set y =   [\cdots]    \in R^n
+          [1]]
+
+Compute b=Ay
+Solve Ax=b
diff --git a/src/approx_derivative.c b/src/approx_derivative.c
new file mode 100644
index 0000000..b33a208
--- /dev/null
+++ b/src/approx_derivative.c
@@ -0,0 +1,38 @@
+#include "lizfcm.h"
+#include <assert.h>
+
+double central_derivative_at(double (*f)(double), double a, double h) {
+  assert(h > 0);
+
+  double x2 = a + h;
+  double x1 = a - h;
+
+  double y2 = (*f)(x2);
+  double y1 = (*f)(x1);
+
+  return (y2 - y1) / (x2 - x1);
+}
+
+double forward_derivative_at(double (*f)(double), double a, double h) {
+  assert(h > 0);
+
+  double x2 = a + h;
+  double x1 = a;
+
+  double y2 = (*f)(x2);
+  double y1 = (*f)(x1);
+
+  return (y2 - y1) / (x2 - x1);
+}
+
+double backward_derivative_at(double (*f)(double), double a, double h) {
+  assert(h > 0);
+
+  double x2 = a;
+  double x1 = a - h;
+
+  double y2 = (*f)(x2);
+  double y1 = (*f)(x1);
+
+  return (y2 - y1) / (x2 - x1);
+}
diff --git a/src/lin.c b/src/lin.c
new file mode 100644
index 0000000..2df6f28
--- /dev/null
+++ b/src/lin.c
@@ -0,0 +1,19 @@
+#include "lizfcm.h"
+#include <assert.h>
+
+Line *least_squares_lin_reg(Array_double *x, Array_double *y) {
+  assert(x->size == y->size);
+
+  uint64_t n = x->size;
+  double sum_x = sum_v(x);
+  double sum_y = sum_v(y);
+  double sum_xy = dot_v(x, y);
+  double sum_xx = dot_v(x, x);
+  double denom = ((n * sum_xx) - (sum_x * sum_x));
+
+  Line *line = malloc(sizeof(Line));
+  line->m = ((sum_xy * n) - (sum_x * sum_y)) / denom;
+  line->a = ((sum_y * sum_xx) - (sum_x * sum_xy)) / denom;
+
+  return line;
+}
diff --git a/src/maceps.c b/src/maceps.c
new file mode 100644
index 0000000..23bc9db
--- /dev/null
+++ b/src/maceps.c
@@ -0,0 +1,28 @@
+#include "lizfcm.h"
+#include <math.h>
+
+float smaceps() {
+  float one = 1.0;
+  float machine_epsilon = 1.0;
+  float one_approx = one + machine_epsilon;
+
+  while (fabsf(one_approx - one) > 0) {
+    machine_epsilon /= 2;
+    one_approx = one + machine_epsilon;
+  }
+
+  return machine_epsilon;
+}
+
+double dmaceps() {
+  double one = 1.0;
+  double machine_epsilon = 1.0;
+  double one_approx = one + machine_epsilon;
+
+  while (fabs(one_approx - one) > 0) {
+    machine_epsilon /= 2;
+    one_approx = one + machine_epsilon;
+  }
+
+  return machine_epsilon;
+}
diff --git a/src/main.c b/src/main.c
new file mode 100644
index 0000000..6bb704a
--- /dev/null
+++ b/src/main.c
@@ -0,0 +1,7 @@
+#include "lizfcm.h"
+#include <stdio.h>
+
+int main() {
+  printf("hello, from lizfcm!\n");
+  return 0;
+}
diff --git a/src/vector,least-squares-reg.lisp b/src/vector,least-squares-reg.lisp
deleted file mode 100644
index 1c7272c..0000000
--- a/src/vector,least-squares-reg.lisp
+++ /dev/null
@@ -1,2 +0,0 @@
-(in-package :lizfcm.vector)
-
diff --git a/src/vector.c b/src/vector.c
new file mode 100644
index 0000000..61692e1
--- /dev/null
+++ b/src/vector.c
@@ -0,0 +1,90 @@
+#include "lizfcm.h"
+#include <assert.h>
+#include <float.h>
+#include <math.h>
+#include <stdio.h>
+
+double l2_norm(Array_double *v) {
+  double norm = 0;
+  for (size_t i = 0; i < v->size; ++i)
+    norm += v->data[i] * v->data[i];
+  return sqrt(norm);
+}
+
+double l1_norm(Array_double *v) {
+  double sum = 0;
+  for (size_t i = 0; i < v->size; ++i)
+    sum += fabs(v->data[i]);
+  return sum;
+}
+
+double linf_norm(Array_double *v) {
+  double max = -DBL_MAX;
+  for (size_t i = 0; i < v->size; ++i)
+    max = c_max(v->data[i], max);
+  return max;
+}
+
+Array_double *minus_v(Array_double *v1, Array_double *v2) {
+  assert(v1->size == v2->size);
+
+  Array_double *sub = InitArrayWithSize(double, v1->size, 0);
+  for (size_t i = 0; i < v1->size; i++)
+    sub->data[i] = v1->data[i] - v2->data[i];
+  return sub;
+}
+
+double sum_v(Array_double *v) {
+  double sum = 0;
+  for (size_t i = 0; i < v->size; i++)
+    sum += v->data[i];
+  return sum;
+}
+
+Array_double *add_v(Array_double *v1, Array_double *v2) {
+  assert(v1->size == v2->size);
+
+  Array_double *sum = InitArrayWithSize(double, v1->size, 0);
+  for (size_t i = 0; i < v1->size; i++)
+    sum->data[i] = v1->data[i] + v2->data[i];
+  return sum;
+}
+
+double dot_v(Array_double *v1, Array_double *v2) {
+  assert(v1->size == v2->size);
+
+  double dot = 0;
+  for (size_t i = 0; i < v1->size; i++)
+    dot += v1->data[i] * v2->data[i];
+  return dot;
+}
+
+double l2_distance(Array_double *v1, Array_double *v2) {
+  Array_double *minus = minus_v(v1, v2);
+  double dist = l2_norm(minus);
+  free(minus);
+  return dist;
+}
+
+double l1_distance(Array_double *v1, Array_double *v2) {
+  Array_double *minus = minus_v(v1, v2);
+  double dist = l1_norm(minus);
+  free(minus);
+  return dist;
+}
+
+double linf_distance(Array_double *v1, Array_double *v2) {
+  Array_double *minus = minus_v(v1, v2);
+  double dist = linf_norm(minus);
+  free(minus);
+  return dist;
+}
+
+void format_vector_into(Array_double *v, char *s) {
+  sprintf(s, "");
+  if (v->size == 0)
+    sprintf(s, "empty");
+
+  for (size_t i = 0; i < v->size; ++i)
+    sprintf(s, "%s %f,", s, v->data[i]);
+}
diff --git a/test/main.c b/test/main.c
new file mode 100644
index 0000000..b0faf97
--- /dev/null
+++ b/test/main.c
@@ -0,0 +1,54 @@
+#include "lizfcm.h"
+#include <stdio.h>
+#include <stdlib.h>
+
+double f(double x) { return (x - 1) / (x + 1); }
+
+int main() {
+  printf("smaceps(): %.10e\n", smaceps());
+  printf("dmaceps(): %.10e\n", dmaceps());
+
+  Array_double *v = InitArray(double, {3, 1, -4, 1, 5, -9, 3});
+  Array_double *w = InitArray(double, {-2, 7, 1, -8, -2, 8, 5});
+
+  char v_s[256];
+  char w_s[256];
+  format_vector_into(v, v_s);
+  format_vector_into(w, w_s);
+
+  printf("v: %s\n", v_s);
+  printf("w: %s\n", w_s);
+  printf("l1_norm(v): %f\n", l1_norm(v));
+  printf("l2_norm(v): %f\n", l2_norm(v));
+  printf("linf_norm(v): %f\n", linf_norm(v));
+
+  printf("l1_dist(v, w): %f\n", l1_distance(v, w));
+  printf("l2_dist(v, w): %f\n", l2_distance(v, w));
+  printf("linf_dist(v, w): %f\n", linf_distance(v, w));
+
+  double h = 0.001;
+  printf("approx f'(1) w/ c.d.: %f\n", central_derivative_at(&f, 1, h));
+  printf("approx f'(1) w/ fw.d.: %f\n", forward_derivative_at(&f, 1, h));
+  printf("approx f'(1) w/ bw.d.: %f\n", backward_derivative_at(&f, 1, h));
+
+  v = InitArray(double, {1, 2, 3, 4, 5});
+  w = InitArray(double, {2, 3, 4, 5, 6});
+  format_vector_into(v, v_s);
+  format_vector_into(w, w_s);
+  printf("v: %s\n", v_s);
+  printf("w: %s\n", w_s);
+
+  Line *line = least_squares_lin_reg(v, w);
+  printf("least_squares_lin_reg(v, w): (%f)x + %f\n", line->m, line->a);
+
+  v = InitArray(double, {1, 2, 3, 4, 5, 6, 7});
+  w = InitArray(double, {0.5, 3, 2, 3.5, 5, 6, 7.5});
+  format_vector_into(v, v_s);
+  format_vector_into(w, w_s);
+  printf("v: %s\n", v_s);
+  printf("w: %s\n", w_s);
+  line = least_squares_lin_reg(v, w);
+  printf("least_squares_lin_reg(v, w): (%f)x + %f\n", line->m, line->a);
+
+  return 0;
+}