实验6——力导向图

#include <bits/stdc++.h>
#include "opencv2/core.hpp"
#include "opencv2/imgproc.hpp"
#include "opencv2/video.hpp"
#include "opencv2/objdetect.hpp"
#include "opencv2/imgcodecs.hpp"
#include "opencv2/highgui.hpp"
#include "opencv2/ml.hpp"
#define inf 2333333333333333
#define N 100010
#define p(a) putchar(a)
#define For(i,a,b) for(int i=a;i<=b;++i)
//by war
//2020.11.29
using namespace std;
using namespace cv;
int n=15,T=100,k,t,d_t=0.1,s;
int mp[100][100];
double L=10.0,k_r=5.0,k_s=6.0,dx,dy,disSqu,dist,force,fx,fy;
Mat image;
double alpha = 0.05;
struct Node{
    double force_x;
    double force_y;
    double x;
    double y;
}a[N];

void in(int &x){
    int y=1;char c=getchar();x=0;
    while(c<'0'||c>'9'){if(c=='-')y=-1;c=getchar();}
    while(c<='9'&&c>='0'){ x=(x<<1)+(x<<3)+c-'0';c=getchar();}
    x*=y;
}
void o(int x){
    if(x<0){p('-');x=-x;}
    if(x>9)o(x/10);
    p(x%10+'0');
}

signed main(){
    srand(time(0));
    For(i,1,n){
        a[i].x=rand()%500+1;
        a[i].y=rand()%500+1;
    }
    For(i,1,n){
        k=rand()%n+1;
        while(k--){
            t=i;
            while(t==i) t=rand()%n+1;
            mp[i][t]=1;
            mp[t][i]=1;
        }
    }
    while(T--){
        For(i,1,n) a[i].force_x=a[i].force_y=0;
        For(i,1,n)
            For(j,i+1,n){
                dx=a[j].x-a[i].x;
                dy=a[j].y-a[i].y;
                if(dx!=0 || dy!=0){
                    disSqu=dx*dx+dy*dy;
                    dist=sqrt(disSqu);
                    force=k_r/disSqu;
                    fx=force*dx/dist;
                    fy=force*dy/dist;
                    a[i].force_x-=fx;
                    a[i].force_y-=fy;
                    a[j].force_x+=fx;
                    a[j].force_y+=fy;
                }
            }
        For(i,1,n)
            For(j,i+1,n){
                if(!mp[i][j]) continue;
                dx=a[j].x-a[i].x;
                dy=a[j].y-a[i].y;
                if(dx!=0 || dy!=0){
                    disSqu=dx*dx+dy*dy;
                    dist=sqrt(disSqu);
                    force=k_s/(dist-L);
                    fx=force*dx/dist;
                    fy=force*dy/dist;
                    a[i].force_x+=fx;
                    a[i].force_y+=fy;
                    a[j].force_x-=fx;
                    a[j].force_y-=fy;
                }
            }
        For(i,1,n){
            dx=d_t*a[i].force_x;
            dy=d_t*a[i].force_y;
            disSqu=dx*dx+dy*dy;
            if(dist>500){
                s=sqrt(500/disSqu);
                dx*=s;
                dy*=s;
            }
            a[i].x+=dx;
            a[i].y+=dy;
        }
    }
    image=Mat(1000, 1000, CV_8UC3, Scalar(0));
    For(i,1,n){
        circle(image, Point(a[i].y,a[i].x), 0.5, Scalar(255, 0, 0), FILLED, 0);
        For(j,i+1,n){
            if(!mp[i][j]) continue;
            line(image, Point(a[i].y,a[i].x), Point(a[j].y,a[j].x), Scalar(0, 255,0), 2, 8, 0);
        }
    }
    imshow("力导向图", image);
    waitKey(0);
    return 0;
}