天天看点

数值计算——施密特正交化

施密特正交化单位正交基C++代码

使用施密特正交化得出单位正交基的完整代码

#include <time.h>
#include <iostream>
#include <string.h>
#include <stdio.h>
#include <vector>
#include <math.h>

using namespace std;
const int N = 5000;
typedef double Type;

Type A[N][N];

void schmidt(Type A[][N], int n, int k, int l); /*施密特正交化*/
void Print(Type L[][N], int n);//输入
void united(Type A[][N], int n);//向量单位化
int Random(int start, int end);/*随机数*/

int Random(int start, int end)
{
	int dis = end - start;
	return rand() % dis + start;
}
void Print(Type L[][N], int n)
{
	for (int i = 0; i < n; i++)
	{
		for (int j = 0; j < n; j++)
			cout << L[i][j] << " ";
		cout << endl;
	}
	cout << endl;
}

void united(Type A[][N], int n) {
	vector<Type> B;
	Type y = 0;
	for (int i = 0; i < n; i++) {
		for (int l=0;l<n;l++)
			y += A[i][l] * A[i][l];
		y = sqrt(y);
		B.push_back(y);
	}
	for (int i = 0; i < n; i++) {
		for (int j = 0; j < n; j++) {
			A[i][j] /= B[i];
		}
	}
	B.clear();
}
void schmidt(Type A[][N],int n,int k,int l) {//施密特正交化
	vector<Type> B;
	double num = 0;
	int t=0;
	//(a.b)/|b|
	for (int i = 0; i < n; i++) {
		t += A[k][i] * A[l][i];
	}
	for (int i = 0; i < n; i++) {
		Type y=0;
		y = t * A[l][i];
		B.push_back(y);
	}
	for (int i = 0; i < n; i++) {
		num += A[l][i] * A[l][i];
	}
	for (int i = 1; i < n ; i++) {
			A[k][i] -= (B[i] / num);
	}
	B.clear();
	if (l+1 < k) {
		schmidt(A, n, k, l + 1);
	}
	else if (k + 1 < n) {
		schmidt(A,  n, k + 1, 0);
	}
}

int main()
{
	srand((int)time(0));
	double duration;
	clock_t start, finish;
	int n;
	cout << "输入矩阵大小";
	cin >> n;
	for (int i = 0; i < n; i++)
	{
		for (int j = 0; j < n; j++) {
			/*if (i == j) {
				A[i][i] = 1;
				continue;
			}
			A[i][j] = 0;*/
			A[i][j] = Random(0, 10);
		}
	}
	//Print(A, n);
	start = clock();//获取当前时间
	schmidt(A, n,1,0);
	//Print(A, n);
	united(A, n);
	//Print(A, n);
	printf("完成!");
	finish = clock();//获取当前时间
	duration = (double)(finish - start);//计算程序运行时间
	cout << "时间为:";
	cout << duration << endl;
	return 0;
}

           

继续阅读