Skip to content

Commit e82ed09

Browse files
committed
Merge branch 'master' into missingexamples
2 parents c7dc32e + a47bdee commit e82ed09

File tree

5 files changed

+127
-28
lines changed

5 files changed

+127
-28
lines changed

chapters/FFT/code/julia/fft.jl

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -5,7 +5,7 @@ function DFT(x)
55
# We want two vectors here for real space (n) and frequency space (k)
66
n = 0:N-1
77
k = n'
8-
transform_matrix = exp.(-2im * pi *n *k / N)
8+
transform_matrix = exp.(-2im*pi*n*k/N)
99
return transform_matrix*x
1010

1111
end
@@ -22,7 +22,6 @@ function cooley_tukey(x)
2222
x_even = x[2]
2323
end
2424
n = 0:N-1
25-
#n = n'
2625
half = div(N,2)
2726
factor = exp.(-2im*pi*n/N)
2827
return vcat(x_odd + x_even .* factor[1:half],
Lines changed: 114 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,114 @@
1+
#include <math.h>
2+
#include <stddef.h>
3+
#include <stdio.h>
4+
#include <stdlib.h>
5+
#include <string.h>
6+
7+
struct point {
8+
double x, y;
9+
};
10+
11+
int cmp_points(const void *a, const void *b) {
12+
struct point* pa = (struct point*) a;
13+
struct point* pb = (struct point*) b;
14+
15+
if (pa->y > pb->y) {
16+
return 1;
17+
} else if (pa->y < pb->y) {
18+
return -1;
19+
} else {
20+
return 0;
21+
}
22+
}
23+
24+
double ccw(struct point a, struct point b, struct point c) {
25+
return (b.x - a.x) * (c.y - a.y) - (b.y - a.y) * (c.x - a.x);
26+
}
27+
28+
double polar_angle(struct point origin, struct point p) {
29+
return atan2(p.y - origin.y, p.x - origin.x);
30+
}
31+
32+
void polar_angles_sort(struct point *points, struct point origin, size_t size) {
33+
if (size < 2) {
34+
return;
35+
}
36+
37+
double pivot_angle = polar_angle(origin, points[size / 2]);
38+
39+
int i = 0;
40+
int j = size - 1;
41+
while (1) {
42+
while (polar_angle(origin, points[i]) < pivot_angle) {
43+
i++;
44+
}
45+
while (polar_angle(origin, points[j]) > pivot_angle) {
46+
j--;
47+
}
48+
49+
if (i >= j) {
50+
break;
51+
}
52+
53+
struct point tmp = points[i];
54+
points[i] = points[j];
55+
points[j] = tmp;
56+
57+
i++;
58+
j--;
59+
}
60+
61+
polar_angles_sort(points, origin, i);
62+
polar_angles_sort(points + i, origin, size - i);
63+
}
64+
65+
size_t graham_scan(struct point *points, size_t size) {
66+
qsort(points, size, sizeof(struct point), cmp_points);
67+
polar_angles_sort(points, points[0], size);
68+
69+
struct point tmp_points[size + 1];
70+
memcpy(tmp_points + 1, points, size * sizeof(struct point));
71+
tmp_points[0] = tmp_points[size];
72+
73+
size_t m = 1;
74+
for (size_t i = 2; i <= size; ++i) {
75+
while (ccw(tmp_points[m - 1], tmp_points[m], tmp_points[i]) <= 0) {
76+
if (m > 1) {
77+
m--;
78+
continue;
79+
} else if (i == size) {
80+
break;
81+
} else {
82+
i++;
83+
}
84+
}
85+
86+
m++;
87+
struct point tmp = tmp_points[i];
88+
tmp_points[i] = tmp_points[m];
89+
tmp_points[m] = tmp;
90+
}
91+
92+
memcpy(points, tmp_points + 1, size * sizeof(struct point));
93+
94+
return m;
95+
}
96+
97+
int main() {
98+
struct point points[] = {{2.0, 1.9}, {1.0, 1.0}, {2.0, 4.0}, {3.0, 1.0},
99+
{2.0, 0.0}};
100+
101+
printf("Points:\n");
102+
for (size_t i = 0; i < 5; ++i) {
103+
printf("(%f,%f)\n", points[i].x, points[i].y);
104+
}
105+
106+
size_t hull_size = graham_scan(points, 5);
107+
108+
printf("\nHull:\n");
109+
for (size_t i = 0; i < hull_size; ++i) {
110+
printf("(%f,%f)\n", points[i].x, points[i].y);
111+
}
112+
113+
return 0;
114+
}

chapters/computational_geometry/gift_wrapping/graham_scan/code/julia/graham.jl

Lines changed: 3 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -7,26 +7,6 @@ function dist(point1::Point, point2::Point)
77
return sqrt((point1.x - point2.x)^2 + (point1.y - point2.y)^2)
88
end
99

10-
function graham_angle(point1::Point, point2::Point, point3::Point)
11-
# Find distances between all points
12-
a = dist(point3, point2)
13-
b = dist(point3, point1)
14-
c = dist(point1, point2)
15-
16-
ret_angle = acos((b*b - a*a - c*c)/(2*a*c))
17-
18-
if(sign(point1.x - point2.x) != sign(point1.x - point3.x))
19-
ret_angle += 0.5*pi
20-
end
21-
22-
if (isnan(ret_angle))
23-
exit(1)
24-
end
25-
26-
return ret_angle
27-
28-
end
29-
3010
function ccw(a::Point, b::Point, c::Point)
3111
return ((b.x - a.x)*(c.y - a.y) - (b.y - a.y)*(c.x - a.x))
3212
end
@@ -38,7 +18,8 @@ function graham_scan(points::Vector{Point})
3818
sort!(points, by = item -> item.y)
3919

4020
# Sort all other points according to angle with that point
41-
other_points = sort(points[2:end], by = item -> graham_angle(Point(points[1].x - 1,points[1].y), points[1], item))
21+
other_points = sort(points[2:end], by = item -> atan2(item.y - points[1].y,
22+
item.x - points[1].x))
4223

4324
# Place points sorted by angle back into points vector
4425
for i in 1:length(other_points)
@@ -69,6 +50,7 @@ function graham_scan(points::Vector{Point})
6950
end
7051

7152
function main()
53+
# This hull is just a simple test so we know what the output should be
7254
points = [Point(2,1.9), Point(1, 1), Point(2, 4), Point(3, 1), Point(2, 0)]
7355
hull = graham_scan(points)
7456
println(hull)

chapters/computational_geometry/gift_wrapping/graham_scan/graham_scan.md

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,8 @@ We can find whether a rotation is counter-clockwise with trigonometric functions
1414
{% method %}
1515
{% sample lang="jl" %}
1616
[import:30-32, lang:"julia"](code/julia/graham.jl)
17+
{% sample lang="c" %}
18+
[import:24-26, lang:"c_cpp"](code/c/graham.c)
1719
{% endmethod %}
1820

1921
If the output of this function is 0, the points are collinear.
@@ -30,6 +32,8 @@ In the end, the code should look something like this:
3032
{% method %}
3133
{% sample lang="jl" %}
3234
[import:34-69, lang:"julia"](code/julia/graham.jl)
35+
{% sample lang="c" %}
36+
[import:65-95, lang:"c_cpp"](code/c/graham.c)
3337
{% endmethod %}
3438

3539
### Bibliography
@@ -41,6 +45,8 @@ In the end, the code should look something like this:
4145
{% method %}
4246
{% sample lang="jl" %}
4347
[import, lang:"julia"](code/julia/graham.jl)
48+
{% sample lang="c" %}
49+
[import, lang:"c_cpp"](code/c/graham.c)
4450
{% endmethod %}
4551

4652
<script>

chapters/matrix_methods/thomas/thomas.md

Lines changed: 3 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,6 @@
1-
##### Dependencies
2-
* [Gaussian Elimination](gaussian_elimination.md)
3-
41
# Thomas Algorithm
52

6-
As alluded to in the Gaussian Elimination Chapter, the Thomas Algorithm (or TDMA -- Tri-Diagonal Matrix Algorithm) allows for programmers to **massively** cut the computational cost of their code from $$\sim O(n^3) \rightarrow \sim O(n)$$! This is done by exploiting a particular case of Gaussian Elimination, particularly the case where our matrix looks like:
3+
As alluded to in the [Gaussian Elimination chapter](../gaussian_elimination/gaussian_elimination.md), the Thomas Algorithm (or TDMA -- Tri-Diagonal Matrix Algorithm) allows for programmers to **massively** cut the computational cost of their code from $$\sim O(n^3) \rightarrow \sim O(n)$$! This is done by exploiting a particular case of Gaussian Elimination, particularly the case where our matrix looks like:
74

85
$$
96
\left[
@@ -19,7 +16,8 @@ $$
1916

2017
By this, I mean that our matrix is *Tri-Diagonal* (excluding the right-hand side of our system of equations, of course!). Now, at first, it might not be obvious how this helps; however, we may divide this array into separate vectors corresponding to $$a$$, $$b$$, $$c$$, and $$d$$ and then solve for $$x$$ with back-substitution, like before.
2118

22-
In particular, we need to find an optimal scale factor for each row and use that. What is the scale factor? Well, it is the diagonal $$-$$ the multiplicative sum of the off-diagonal elements. In the end, we will update $$c$$ and $$d$$ to be$$c'$$ and $$d'$$ like so:
19+
In particular, we need to find an optimal scale factor for each row and use that. What is the scale factor? Well, it is the diagonal $$-$$ the multiplicative sum of the off-diagonal elements.
20+
In the end, we will update $$c$$ and $$d$$ to be $$c'$$ and $$d'$$ like so:
2321

2422
$$
2523
\begin{align}

0 commit comments

Comments
 (0)