假设我在一个球体上有10个点(随机分布),我想旋转整个系统以确保一个点位于北极。我将如何使用C++来做到这一点?旋转点的球形系统
我去它通过看3D旋转矩阵:
http://en.wikipedia.org/wiki/Rotation_matrix
我转动我的绕X轴点,直到y分量是零,然后绕Y轴直到x分量是零。这应该将问题置于北极或南极吗?
我的代码看起来像这样:
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <iostream>
#include <iomanip>
#include <fstream>
#include <time.h>
#include <stdlib.h>
#include <sstream>
using namespace std;
#define PI 3.14159265358979323846
int main()
{
int a,b,c,f,i,j,k,m,n,s;
double p,Time,Averagetime,Energy,energy,Distance,Length,DotProdForce,Forcemagnitude,
ForceMagnitude[101],Force[101][4],E[1000001],En[501],x[101][4],y[101][4],
Dist[101][101],Epsilon,z[101],theta,phi;
/* set the number of points */
n=10;
/* check that there are no more than 100 points */
if(n>100){
cout << n << " is too many points for me :-(\n";
exit(0);
}
/* reset the random number generator */
srand((unsigned)time(0));
for (i=1;i<=n;i++){
x[i][1]=((rand()*1.0)/(1.0*RAND_MAX)-0.5)*2.0;
x[i][2]=((rand()*1.0)/(1.0*RAND_MAX)-0.5)*2.0;
x[i][3]=((rand()*1.0)/(1.0*RAND_MAX)-0.5)*2.0;
Length=sqrt(pow(x[i][1],2)+pow(x[i][2],2)+pow(x[i][3],2));
for (k=1;k<=3;k++){
x[i][k]=x[i][k]/Length;
}
}
/* calculate the energy */
Energy=0.0;
for(i=1;i<=n;i++){
for(j=i+1;j<=n;j++){
Distance=sqrt(pow(x[i][1]-x[j][1],2)+pow(x[i][2]-x[j][2],2)
+pow(x[i][3]-x[j][3],2));
Energy=Energy+1.0/Distance;
}
}
cout << fixed << setprecision(10) << "energy=" << Energy << "\n";
/* Save Values so far */
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
y[i][j]=x[i][j];
}
}
/* Choose each point in turn and make it the north pole note this is what the while loop os for, but have set it to a<2 to just look at first point */
a=1;
b=0;
c=0;
while(a<2){
/* Find theta and phi to rotate points by */
cout << fixed << setprecision(5) << "x[" << a << "][1]=" << x[a][1] <<
" x[" << a << "][2]=" << x[a][2] << " x[" << a << "][3]=" << x[a][3] << "\n";
theta=x[a][2]/x[a][3];
theta=b*PI+atan(theta);
/* Rotate Points by theta around x axis and then by phi around y axis */
for(i=1;i<=n;i++){
x[i][1]=x[i][1];
x[i][2]=x[i][2]*cos(theta)-x[i][3]*sin(theta);
x[i][3]=x[i][2]*sin(theta)+x[i][3]*cos(theta);
}
phi=x[a][1]/x[a][3];
phi=c*PI+atan(phi);
for(i=1;i<=n;i++){
x[i][1]=x[i][1]*cos(phi)-x[i][3]*sin(phi);
x[i][2]=x[i][2];
x[i][3]=x[i][1]*sin(phi)+x[i][3]*cos(phi);
}
cout << fixed << setprecision(5) << "x[" << a << "][1]=" << x[a][1] <<
" x[" << a << "][2]=" << x[a][2] << " x[" << a << "][3]=" << x[a][3] << "\n";
if(x[a][3]==1.0 && x[a][2]==x[a][3]==0)
a=a+1;
else if(b==0 && c==0)
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
x[i][j]=y[i][j];
c=1;
}
}
else if(b==0 && c==1)
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
x[i][j]=y[i][j];
b=1;
c=0;
}
}
else if(b==1 && c==0)
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
x[i][j]=y[i][j];
c=1;
}
}
else if(b==1 && c==1)
break;
}
energy=0.0;
for(i=1;i<=n;i++){
for(j=i+1;j<=n;j++){
Distance=sqrt(pow(x[i][1]-x[j][1],2)+pow(x[i][2]-x[j][2],2)
+pow(x[i][3]-x[j][3],2));
energy=energy+1.0/Distance;
}
}
cout << fixed << setprecision(10) << "ENERGY=" << energy << "\n";
cout << fixed << setprecision(5) << "x[" << a << "][1]=" << x[a][1] <<
" x[" << a << "][2]=" << x[a][2] << " x[" << a << "][3]=" << x[a][3] << "\n";
/* Output to help with gnuin.txt */
ofstream File4 ("mypoints");
for(i=1;i<=n;i++){
File4 << x[i][1] << " " << x[i][2] << " " << x[i][3] << "\n";
}
File4.close();
return 0;
}
好了,这里有问题的负载,就像if语句(线103)真的不应该有平等的条件双重因为它会从来没有工作,但我可以稍后使用间接比较(一些epsilon的东西)。我真正的疑问是为什么轮转,即使它在所有点上作用,都会将球体上的点取下来? (正如你所看到的,这些值已经过标准化,使它们全部在第38行的单位球面上)。
注意:b,c的东西是检查点是在北极还是南极。
在这一点上,你的问题实际上是低估的:自由度不受约束导致整个类的最终状态。这可能或可能不是一个问题。 – dmckee
您必须添加一些约束条件。只需选取一个随机点并沿着与该点相交的平面旋转整个系统点,即球体的中心和球体顶部的虚拟点。 – Blender