pzr1988 于 2013.07.08 21:10 提问

#include
#define DIST(X,Y,A,B) DA=(X)-(A);DB=(Y)-(B);D=std::sqrt(DA*DA+DB*DB);C=std::max(1.0-(D/H)*(D/H)*(D/H),0.0)*100.0f;
int NRX=ceil(X/H),NRY=ceil(Y/H);
std::vector FIRST;
std::vectorstd::pair<std::pair<std::pair<double,double,std::pair >,int> > PARTICLE;
void step(double dt,int i){
FIRST.assign(NRX*NRY,-1);
for(size_t i=0;i PARTICLE[i].second=-1;
double D,DA,DB,C;
DIST(PARTICLE[i].first.first.first,PARTICLE[i].first.first.second,X/2.0,Y/2.0)
int& f=FIRST[std::floor(PARTICLE[i].first.first.first/H)*NRY+std::floor(PARTICLE[i].first.first.second/H)];
PARTICLE[i].second=f;
f=i;
}
for(size_t i=0;i PARTICLE[i].first.second.second+=G*dt;
int HX=std::floor(PARTICLE[i].first.first.first/H),HY=std::floor(PARTICLE[i].first.first.second/H);
for(int XX=std::max(HX-1,0);XX for(int YY=std::max(HY-1,0);YY for(int OFF=FIRST[XX*NRY+YY];OFF!=-1;OFF=PARTICLE[OFF].second){
double D,DA,DB,C;
DIST(PARTICLE[i].first.first.first,PARTICLE[i].first.first.second,PARTICLE[OFF].first.first.first,PARTICLE[OFF].first.first.second)
if(OFF != i && C > 0.0){PARTICLE[i].first.second.first+=DA*C*dt/D;PARTICLE[i].first.second.second+=DB*C*dt/D;}
}
}
for(size_t i=0;i PARTICLE[i].first.second.first*=D;
PARTICLE[i].first.second.second*=D;
PARTICLE[i].first.first.first+=PARTICLE[i].first.second.first*dt;
PARTICLE[i].first.first.second+=PARTICLE[i].first.second.second*dt;
}
}
void init(){
for(double XX=0;XX for(double YY=0;YY double D,DA,DB,C;
DIST(XX,YY,X/2.0,Y/2.0);
if(D PARTICLE.push_back(std::make_pair(std::make_pair(std::make_pair(XX,YY),std::make_pair(0.0,0.0)),-1));
}
}
}
#include
int main()
{
init();
HWND myconsole = GetConsoleWindow();
HDC mydc = GetDC(myconsole);
SelectObject(mydc,CreateSolidBrush(RGB(0,0,0)));
for(int i=0;i<100000;i++){
step(0.002f,i);
if(i%100==0){
Rectangle(mydc,0,0,X*10,Y*10);
for(size_t i=0;i<PARTICLE.size();i++)SetPixel(mydc,PARTICLE[i].first.first.first*10,10*Y-PARTICLE[i].first.first.second*10,RGB(255,255,255));
}
}
ReleaseDC(myconsole, mydc);
system("pause");
}

2个回答

pzr1988   2013.07.09 09:37
``````#include <vector>
#define DIST(X,Y,A,B) DA=(X)-(A);DB=(Y)-(B);D=std::sqrt(DA*DA+DB*DB);C=std::max(1.0-(D/H)*(D/H)*(D/H),0.0)*100.0f;
int NRX=ceil(X/H),NRY=ceil(Y/H);
std::vector<int> FIRST;
std::vector<std::pair<std::pair<std::pair<double,double>,std::pair<double,double> >,int> > PARTICLE;
void step(double dt,int i){
FIRST.assign(NRX*NRY,-1);
#pragma omp parallel for
for(size_t i=0;i<PARTICLE.size();i++){
PARTICLE[i].second=-1;
double D,DA,DB,C;
DIST(PARTICLE[i].first.first.first,PARTICLE[i].first.first.second,X/2.0,Y/2.0)
int& f=FIRST[std::floor(PARTICLE[i].first.first.first/H)*NRY+std::floor(PARTICLE[i].first.first.second/H)];
PARTICLE[i].second=f;
f=i;
}
#pragma omp parallel for
for(size_t i=0;i<PARTICLE.size();i++){
PARTICLE[i].first.second.second+=G*dt;
int HX=std::floor(PARTICLE[i].first.first.first/H),HY=std::floor(PARTICLE[i].first.first.second/H);
for(int XX=std::max(HX-1,0);XX<=std::min(HX+1,NRX-1);XX++)
for(int YY=std::max(HY-1,0);YY<=std::min(HY+1,NRY-1);YY++)
for(int OFF=FIRST[XX*NRY+YY];OFF!=-1;OFF=PARTICLE[OFF].second){
double D,DA,DB,C;
DIST(PARTICLE[i].first.first.first,PARTICLE[i].first.first.second,PARTICLE[OFF].first.first.first,PARTICLE[OFF].first.first.second)
if(OFF != i && C > 0.0){PARTICLE[i].first.second.first+=DA*C*dt/D;PARTICLE[i].first.second.second+=DB*C*dt/D;}
}
}
#pragma omp parallel for
for(size_t i=0;i<PARTICLE.size();i++){
PARTICLE[i].first.second.first*=D;
PARTICLE[i].first.second.second*=D;
PARTICLE[i].first.first.first+=PARTICLE[i].first.second.first*dt;
PARTICLE[i].first.first.second+=PARTICLE[i].first.second.second*dt;
}
}
void init(){
for(double XX=0;XX<X;XX+=H*0.5f)
for(double YY=0;YY<Y;YY+=H*0.5f){
double D,DA,DB,C;
DIST(XX,YY,X/2.0,Y/2.0);
PARTICLE.push_back(std::make_pair(std::make_pair(std::make_pair(XX,YY),std::make_pair(0.0,0.0)),-1));
}
}
}
#include <Windows.h>
int main()
{
init();
HWND myconsole = GetConsoleWindow();
HDC mydc = GetDC(myconsole);
SelectObject(mydc,CreateSolidBrush(RGB(0,0,0)));
for(int i=0;i<100000;i++){
step(0.002f,i);
if(i%100==0){
Rectangle(mydc,0,0,X*10,Y*10);
for(size_t i=0;i<PARTICLE.size();i++)SetPixel(mydc,PARTICLE[i].first.first.first*10,10*Y-PARTICLE[i].first.first.second*10,RGB(255,255,255));
}
}
ReleaseDC(myconsole, mydc);
system("pause");
}
``````
John_ToString   2015.12.05 16:35