Skip to content

Instantly share code, notes, and snippets.

@miettal
Created May 13, 2014 15:32
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save miettal/8a60decbdd2936df1746 to your computer and use it in GitHub Desktop.
Save miettal/8a60decbdd2936df1746 to your computer and use it in GitHub Desktop.
グラムシュミットの正規直交化プログラムらしい.授業でやった.内容は忘れた.
#include <stdio.h>
#include <math.h>
#define MAX_DIM 128
#define DIM 3
typedef double vector[MAX_DIM];
double vectorInnterProduct(int dim, vector vec1, vector vec2)
{
int i;
double ret = 0;
for (i = 0; i < dim; i++) {
ret += vec1[i] * vec2[i];
}
return ret;
}
void vectorMul(int dim, vector vec1, double d, vector vec2)
{
int i;
for (i = 0; i < dim; i++) {
vec2[i] = d * vec1[i];
}
}
void vectorSub(int dim, vector vec1, vector vec2, vector vec3)
{
int i;
for (i = 0; i < dim; i++) {
vec3[i] = vec1[i] - vec2[i];
}
}
void vectorAdd(int dim, vector vec1, vector vec2, vector vec3)
{
int i;
for (i = 0; i < dim; i++) {
vec3[i] = vec1[i] + vec2[i];
}
}
double vectorNorm(int dim, vector vec)
{
return pow(vectorInnterProduct(dim, vec, vec), 1.0/2.0);
}
void vectorPrint(int dim, vector vec)
{
int i;
for (i = 0; i < dim; i++) {
printf("%f ", vec[i]);
}
puts("");
return;
}
void vectorOrthonormalize(int dim, vector v[], vector w[])
{
int n, i;
vectorMul(dim, v[0], 1.0/vectorNorm(dim, v[0]), w[0]);
for (n = 1; n < dim; n++) {
vector t;
vector s = {0};
vector w_;
for (i = 0; i < n; i++) {
vectorMul(dim, w[i], vectorInnterProduct(dim, w[i], v[n]), t);
vectorAdd(dim, t, s, s);
}
vectorSub(dim, v[n], s, w_);
vectorMul(dim, w_, 1.0/vectorNorm(dim, w_), w[n]);
}
return;
}
int main(int argc, const char *argv[])
{
int i, j;
int dim;
vector v[MAX_DIM];
vector w[MAX_DIM];
printf("次元数:"); scanf("%d", &dim);
for (i = 0; i < dim; i++) {
for (j = 0; j < dim; j++) {
printf("v[%d][%d]:", i, j); scanf("%lf", &v[j][i]);
}
}
vectorOrthonormalize(dim, v, w);
for (i = 0; i < DIM; i++) {
vectorPrint(dim, w[i]);
}
return 0;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment