diff --git a/graph.py b/graph.py
index cafe73cc418764c315b4e9576ab2fcd1d1c9dd02..1b312b3d4b9dd6df23ad72b0d07473b886f8617c 100644
--- a/graph.py
+++ b/graph.py
@@ -2,29 +2,34 @@ import sys
import matplotlib.pyplot as plt
import numpy as np
-x = []
-y = []
-a = 0
-b = 0
-
-with open('points.txt') as f:
- lines = f.readlines()
-
- z = lines[0].split(",")
- a = float(z[0])
- b = float(z[1])
-
- for line in lines[1:]:
+x = [[], [], []]
+y = [[], [], []]
+a = []
+b = []
+
+fig, axes = plt.subplots(3, 1, figsize=(18, 15))
+plt.get_current_fig_manager().set_window_title('Graphics')
+
+for i in range(3):
+ with open('G' + str(i+1) + '.txt') as f:
+ for line in f.readlines():
+ z = line.split(",")
+ x[i].append(float(z[0]))
+ y[i].append(float(z[1]))
+
+with open('ab.txt') as f:
+ for line in f.readlines():
z = line.split(",")
- x.append(float(z[0]))
- y.append(float(z[1]))
+ a.append(float(z[0]))
+ b.append(float(z[1]))
-fig, ax = plt.subplots()
+titles = [(2, r'$G1 \bigcup G2 \longmapsto G3$'), (0, r'$G2 \bigcup G3 \longmapsto G1$'), (1, r'$G1 \bigcup G3 \longmapsto G2$')]
# 1st point (0, b) and 2nd point (b, a+b)
-xLine, yLine = [0, 1], [b, a+b]
-plt.plot(xLine, yLine, marker = 'o')
-
-ax.plot(x, y, 'o', color='tab:brown')
+for i in range(3):
+ axes[i].set_title(titles[i][1])
+ axes[i].axline([0, b[i]], [1, a[i] + b[i]])
+ axes[i].plot(x[titles[i][0]], y[titles[i][0]], '.', color='tab:green')
+ axes[i].grid(True)
plt.show()
\ No newline at end of file
diff --git a/main.c b/main.c
index 03adc33d14f20bd06b1033007723c6313502f755..80df2b973f58beb0a0eede4c27818bb1c5a0b3f2 100644
--- a/main.c
+++ b/main.c
@@ -2,18 +2,21 @@
#include <stdlib.h>
#include <math.h>
#include <time.h>
-//#include <system.h>
#include "vec.h"
#include "opti_num.h"
+#define NB_G 3
+
double randNum();
-void write_to_file(vec ab, vec *points, int N, char *filename);
+int *rand_order_indexes(int N);
+vec **cut_table(int *rand_indexes, int N, vec *points);
+void write_to_file(vec *points, int N, char *filename);
int main() {
srand(time(0));
- int N = 100;
+ int N = NB_G * 200;
vec *points = malloc(sizeof(vec) * N);
double c = randNum(), d = randNum();
@@ -23,26 +26,57 @@ int main() {
vec tmp;
tmp.x = randNum();
- r = randNum() / 100.0;
+ r = randNum() / 10.0;
tmp.y = c * tmp.x + d + r;
points[i] = tmp;
}
- vec ab;
- ab.x = rand() % 10;
- ab.y = rand() % 10;
+ int *rnd_indexes = rand_order_indexes(N);
+
+ vec **groups = cut_table(rnd_indexes, N, points);
- PrintVec(ab);
+ vec *ab = malloc(sizeof(vec) * NB_G);
+ for (int i = 0; i < NB_G; i++) {
+ ab[i].x = randNum();
+ ab[i].y = randNum();
- minimiser(points, N, 0.001, &ab);
- printf("%g, %g\n", c, d);
+ int other = i + 1, to_test = i - 1;
+ if (other >= NB_G) {
+ other = 0;
+ } else if (to_test < 0) {
+ to_test = NB_G - 1;
+ }
- write_to_file(ab, points, N, "points.txt");
+ minimiser(groups[i], groups[other], N / NB_G, 1e-6, &ab[i]);
+ printf("model: a: %g, b: %g\n", ab[i].x, ab[i].y);
+
+ vec test_ab;
+ calc_ab(groups[to_test], NB_G, &test_ab);
+ printf("test: a: %g, b: %g\n\n", test_ab.x, test_ab.y);
+ }
+
+ //minimiser(points, N, 0.001, &ab);
+ printf("c: %g, d: %g\n", c, d);
+
+ for (int i = 0; i < NB_G; i++) {
+ char fname[11];
+ sprintf(fname, "G%d.txt", i + 1);
+ write_to_file(groups[i], N / NB_G, fname);
+ }
+
+ write_to_file(ab, NB_G, "ab.txt");
system("python3 graph.py");
+ free(ab);
free(points);
+ free(rnd_indexes);
+
+ for (int i = 0; i < NB_G; i++) {
+ free(groups[i]);
+ }
+ free(groups);
return EXIT_SUCCESS;
}
@@ -51,15 +85,43 @@ double randNum() {
return rand() / (double)RAND_MAX;
}
-void write_to_file(vec ab, vec *points, int N, char *filename) {
+vec **cut_table(int *rand_indexes, int N, vec *points) {
+ vec **groups = malloc(sizeof(vec*) * NB_G);
+
+ int frac = N / NB_G;
+
+ for (int i = 0; i < NB_G; i++) {
+ groups[i] = malloc(sizeof(vec) * frac);
+ }
+
+ for (int i = 0; i < N; i++) {
+ groups[i / frac][i % frac] = points[rand_indexes[i]];
+ }
+
+ return groups;
+}
+
+int *rand_order_indexes(int N) {
+ int *indexes = calloc(0, sizeof(int) * N);
+
+ for (int i = 1; i < N; i++) {
+ int rnd_index = rand() % N;
+ while (indexes[rnd_index] != 0) {
+ rnd_index = rand() % N;
+ }
+ indexes[rnd_index] = i;
+ }
+
+ return indexes;
+}
+
+void write_to_file(vec *points, int N, char *filename) {
FILE *fp = fopen(filename, "w");
if (fp == NULL) {
return;
}
- fprintf(fp, "%f,%f\n", ab.x, ab.y);
-
for (int i = 0; i < N; i++) {
fprintf(fp, "%f,%f\n", points[i].x, points[i].y);
}
diff --git a/opti_num.c b/opti_num.c
index 91ea76d98e3c8070dc17816165feb79f82c72be3..a263a32761ad92cc4e54f8e6c9b895a0214ef3ba 100644
--- a/opti_num.c
+++ b/opti_num.c
@@ -3,31 +3,48 @@
#include <math.h>
#include <stdio.h>
-void minimiser(vec *vecs, int N, double lambda, vec *ab) {
+void minimiser(vec *vecs1, vec *vecs2, int N, double lambda, vec *ab) {
double sum_x = 0, sum_xx = 0, sum_xy = 0, sum_y = 0;
- for (int i = 0; i < N; i++) {
- sum_x += vecs[i].x;
- sum_xx += vecs[i].x * vecs[i].x;
- sum_xy += vecs[i].x * vecs[i].y;
- sum_y += vecs[i].y;
+ for (int i = 0; i < 2 * N; i++) {
+ if (i < N) {
+ sum_x += vecs1[i].x;
+ sum_xx += vecs1[i].x * vecs1[i].x;
+ sum_xy += vecs1[i].x * vecs1[i].y;
+ sum_y += vecs1[i].y;
+ } else {
+ int tmp_i = i % N;
+
+ sum_x += vecs2[tmp_i].x;
+ sum_xx += vecs2[tmp_i].x * vecs2[tmp_i].x;
+ sum_xy += vecs2[tmp_i].x * vecs2[tmp_i].y;
+ sum_y += vecs2[tmp_i].y;
+ }
}
-
+
vec next_ab = *ab, err_quad;
do {
*ab = next_ab;
err_quad.x = ab->x * sum_xx + ab->y * sum_x - sum_xy;
- err_quad.y = ab->x * sum_x + N * ab->y - sum_y;
+ err_quad.y = ab->x * sum_x + (2 * N) * ab->y - sum_y;
next_ab = Subtract(*ab, Multiply(err_quad, lambda));
- printf("%g, %g\n", next_ab.x, next_ab.y);
} while (GetNorme(Subtract(next_ab, *ab)) > EPSILON);
}
-void calc_ab(double sum_x, double sum_xx, double sum_xy, double sum_y, int N, vec *ab) {
+void calc_ab(vec *points, int N, vec *ab) {
+ double sum_x = 0, sum_xx = 0, sum_xy = 0, sum_y = 0;
+
+ for (int i = 0; i < N; i++) {
+ sum_x += points[i].x;
+ sum_xx += points[i].x * points[i].x;
+ sum_xy += points[i].x * points[i].y;
+ sum_y += points[i].y;
+ }
+
ab->x = (sum_y * sum_x - N * sum_xy) / (sum_x * sum_x - N * sum_xx);
ab->y = (sum_y - ab->x * sum_x) / N;
diff --git a/opti_num.h b/opti_num.h
index 29826ff195861ac3be9f10bc105b635614d4aa27..2524e6b99afe0c29c270dbafc22bea0490fa83c5 100644
--- a/opti_num.h
+++ b/opti_num.h
@@ -4,7 +4,7 @@
#define _OPTI_NUM_H_
#define EPSILON 1e-6
-void minimiser(vec *vecs, int N, double lambda, vec *ab);
-void calc_ab(double sum_x, double sum_xx, double sum_xy, double sum_y, int N, vec *ab);
+void minimiser(vec *vecs1, vec *vecs2, int N, double lambda, vec *ab);
+void calc_ab(vec *points, int N, vec *ab);
#endif