`
暴风雪
  • 浏览: 388107 次
  • 性别: Icon_minigender_2
  • 来自: 杭州
社区版块
存档分类
最新评论

[opengl]弹簧质点模型的opengl实现

 
阅读更多

实现了一种基于弹簧质点模型的布料仿真程序

论文见:http://wenku.baidu.com/link?url=8rBDraTsWTtTL8cFU0vjhXMyv4RF0npjeizz2CQQI4DvvTrxsN3bNOK91_1jRw7TuVadHuds5VnWzQ8CxP0QSOcY6sEUMKUib44crkbil0K

 

对于超弹问题,依然没能很好地解决,希望大神给予指导

#ifndef GLUT_DISABLE_ATEXIT_HACK
#define GLUT_DISABLE_ATEXIT_HACK
#endif
#define GLEW_STATIC
#include <GL/glew.h>
#include <GL/wglew.h>
#include <GL/freeglut.h>
#include <vector>
#include<cstring>
#include <glm/glm.hpp>
#include <glm/gtc/matrix_transform.hpp> //for matrices
#include <glm/gtc/type_ptr.hpp>
//undefine if u want to use the default bending constraint of pbd
#include<iostream>
using namespace std;
//using namespace glm;

GLfloat rtx = 0.0f, rty = 0.0f, rtz = 0.0;
GLfloat sim_step = 0.1;

glm::vec3 Up=glm::vec3(0,1,0), Right, viewDir;
const int GRID_SIZE=10;  //地板砖边长
GLdouble MV[16];
GLint viewport[4];
GLdouble PP[16];
bool isfix[1000];
int dist = -23;

//椭球参数
int iStacks = 30;
int iSlices = 30;
float fRadius = 1;
//float exp1 = 1e-3;

const int numX = 20, numY=20; //一行有numx+1个点
const int total_points = (numX+1)*(numY+1); //总点数
//布料顶点位置 速度
const int fck = 10000;
glm::vec3 Pos[total_points];
glm::vec3 Veloc[total_points];
glm::vec3 force[total_points];
int size = 4;
float hsize = size/2.0f;
float currentTime,newTime,startTime;

const float mass = 1.0/total_points;
const glm::vec3 gvat = glm::vec3(0,-9.8,0);  //重力加速度
const float Spring_K = 2.5;  //弹性系数
const float len0 = 4.0/numX; //单边长度
const float tolera = 1.08;   //弹性限度

void initGL()
{
    //初始化顶点位置
	memset(Pos,0,sizeof(Pos));
	memset(Veloc,0,sizeof(Veloc));
	memset(force,0,sizeof(force));
    //fill in positions
    int count1 = 0;
    int u = numX + 1;
    int v = numY + 1;
	for(int j=0;j<=numY;j++) {
		for(int i=0;i<=numX;i++) {
            Pos[count1++] = glm::vec3( ((float(i)/(u-1)) *2-1)* hsize, size+1, ((float(j)/(v-1) )* size));
			printf("(%.1lf ,%.1lf)",((float(i)/(u-1)) *2-1)* hsize,((float(j)/(v-1) )* size));
		}printf("\n");
		//悬挂点为X[0] 和 X[numX]
	}
	memset(isfix,0,sizeof(isfix));
	isfix[0] = isfix[numX] = 1;

	//处理时间
	startTime = (float)glutGet(GLUT_ELAPSED_TIME);
	cout<<startTime<<"ddddddddddd\n";
	currentTime = startTime;
}

void DrawGrid()
{
	glBegin(GL_LINES);
	glColor3f(0.5f, 0.5f, 0.5f);
	for(int i=-GRID_SIZE;i<=GRID_SIZE;i++)
	{
		glVertex3f((float)i,0,(float)-GRID_SIZE);
		glVertex3f((float)i,0,(float)GRID_SIZE);

		glVertex3f((float)-GRID_SIZE,0,(float)i);
		glVertex3f((float)GRID_SIZE,0,(float)i);
	}
	glEnd();
}

void drawTextile(){
	//draw polygons
	int k = 0;
    glBegin(GL_LINES);
    glColor3f(1,1,1);
	for(int i = 0;i <= numX;i ++){
        for(int j = 0;j <= numY;j ++){
//                cout<<i<<" "<<j<<";";
            if(j != numX){
                glm::vec3 p1 = Pos[k];
                glm::vec3 p2 = Pos[k+1];
                glVertex3f(p1.x,p1.y,p1.z);
                glVertex3f(p2.x,p2.y,p2.z);
//                printf("add %d %d\n",k,k+1);
            }
            if(i != numY){
                glm::vec3 p1 = Pos[k];
                glm::vec3 p2 = Pos[k+numX+1];
                glVertex3f(p1.x,p1.y,p1.z);
                glVertex3f(p2.x,p2.y,p2.z);
            }
            k++;
        }
//        cout<<endl;
	}
    glEnd();

    glPointSize(3);
    glBegin(GL_POINTS);
    glColor3f(1,0,0);
	for(size_t i = 0;i < total_points;i ++){
        glm::vec3 p = Pos[i];
		glVertex3f(p.x,p.y,p.z);
	}
    glEnd();
}

void OnRender(void)
{
    glClear(GL_COLOR_BUFFER_BIT| GL_DEPTH_BUFFER_BIT);
    glLoadIdentity();
    //设置视角
	glTranslatef(0,0,dist);
	glRotatef(15,1,0,0);
	//做出地面
    DrawGrid();
    drawTextile();
//    drawblock();
    glutSwapBuffers();
}
//求出b点对a点的作用力
glm::vec3 SpringForce(int a, int b, int is){
    float inilen;
    if(is == 1)inilen = len0;
    if(is == 2) inilen = len0 * 2;
    if(is == 3) inilen = len0 * 1.414213;
    glm::vec3 res = glm::vec3(0);
    glm::vec3 tmp = Pos[b] - Pos[a];
    float dis = glm::length(tmp);
    res = tmp/dis;
    res *= (dis - inilen);
    res *= Spring_K;
    return res;
}
void ComputeForces(){
    //重力
	for(size_t i=0;i<total_points;i++) {
		force[i] = glm::vec3(0);
		if(!isfix[i])
            force[i] += mass * gvat;
	}
	int i,j,k = 0;
    //结构弹簧 上下左右
    for(i = 0;i <= numY;i++){
        for(j = 0;j <= numX;j ++){
            if(i!=0){//上
                force[k] += SpringForce(k, k - numX - 1, 1);
            }
            if(j!=numX){//右
                force[k] += SpringForce(k, k+1, 1);
            }
            if(i!=numY){//下
                force[k] += SpringForce(k, k + numX + 1, 1);
            }
            if(j!=0){//左
                force[k] += SpringForce(k, k - 1, 1);
            }
            k++;
        }
    }
    //柔性弹簧
    k = 0;
    for(i = 0;i <= numY;i++){
        for(j = 0;j <= numX;j ++){
            if(i>1){//上
                force[k] += SpringForce(k, k - 2*numX - 2, 2);
            }
            if(j<numX-1){//右
                force[k] += SpringForce(k, k+2, 2);
            }
            if(i<numY-1){//下
                force[k] += SpringForce(k, k + 2*numX + 2, 2);
            }
            if(j>1){//左
                force[k] += SpringForce(k, k - 2, 2);
            }
            k++;
        }
    }
    //剪切弹簧
    k = 0;
    for(i = 0;i <= numY;i ++){
        for(j = 0;j <= numX;j ++){
            if(i>0&&j>0){  //左上
                force[k] += SpringForce(k, k - numX - 2, 3);
            }
            if(i>0&&j<numX){  //右上
                force[k] += SpringForce(k, k - numX, 3);
            }
            if(i<numY&&j<numX){  //右下
                force[k] += SpringForce(k, k + numX + 2, 3);
            }
            if(i<numY&&j>0){   //坐下
                force[k] += SpringForce(k, k + numX, 3);
            }
            k ++;
        }
    }
}

void CalcPos(){  //计算新的位置
    newTime = (float) glutGet(GLUT_ELAPSED_TIME);
	float frameTime = newTime-currentTime;
	frameTime/=1000;
	currentTime = newTime;
	glm::vec3 acc = glm::vec3(0);
    //得到的frame就是间隔时间 毫秒
    for(size_t i = 0;i <= total_points;i ++){
        if(isfix[i])continue;
        acc = force[i] / mass;  //得到加速度向量
        Veloc[i] *= 0.98;
        Veloc[i] = Veloc[i] + acc*frameTime;  //得到新的速度值
        Pos[i] = Pos[i] + Veloc[i] * frameTime;
    }
}

void gao(int a,int b,int is){  //直线弹簧 要考虑是否是isfix!
    glm::vec3 tmp = Pos[b] - Pos[a];
    float dis = glm::length(tmp);
    float inilen;
    if(is == 1)inilen = len0;
    if(is == 2) inilen = len0 * 2;
    if(is == 3) inilen = len0 * 1.414213;
    inilen *= tolera;
    if(dis<=inilen)return;
    glm::vec3 tmp1 = tmp*(inilen/dis);
    if(isfix[a]&&isfix[b]){
        return;
    }else if(isfix[a]){
        Pos[b] = Pos[a] + tmp1;
    }else if(isfix[b]){
        Pos[a] = Pos[b] - tmp1;
    }else{
        glm::vec3 mid = (Pos[a] + Pos[b]) * 0.5f;
        Pos[a] = mid - tmp1 * 0.5f;
        Pos[b] = mid + tmp1 * 0.5f;
    }
}

//限制每一根弹簧长度
void DynamicConstrain(){
	int i,j,k = 0;
    for(i = 0;i <= numY;i++){
        for(j = 0;j <= numX;j ++){
            if(j!=numX){//右
                gao(k,k+1,1);
            }
            if(i!=numY){//下
                gao(k,k + numX + 1,1);
            }


            if(j<numX-1){//右二
                gao(k,k + 2, 2);
            }
            if(i<numY-1){//下二
                gao(k, k + 2*numX + 2, 2);
            }


            if(i<numY&&j<numX){  //右下
                gao(k,  k + numX + 2, 3);
            }
            if(i<numY&&j>0){   //左下
                gao(k, k + numX, 3);
            }
            k++;
        }
    }
}
void StepPhysics(){
    ComputeForces();
    CalcPos();
    glutPostRedisplay();
//    DynamicConstrain();
    Sleep(5);
}


void OnReshape(int nw, int nh) {
	glViewport(0,0,nw, nh);
	glMatrixMode(GL_PROJECTION);
	glLoadIdentity();
	gluPerspective(60, (GLfloat)nw / (GLfloat)nh, 1.f, 100.0f);

	glGetIntegerv(GL_VIEWPORT, viewport);
	glGetDoublev(GL_PROJECTION_MATRIX, PP);

	glMatrixMode(GL_MODELVIEW);
}
int main(int argc, char * argv[])
{
	glutInit(&argc, argv);
	glutInitDisplayMode(GLUT_DOUBLE | GLUT_RGBA | GLUT_DEPTH);
    glutInitWindowSize(1024, 1024);
    glutCreateWindow("Mass-Spring Model");

    initGL();
    glutDisplayFunc(OnRender);
    //指定窗口形状变化时的回调函数
	glutReshapeFunc(OnReshape);
	//指定程序空闲时调用函数
    glutIdleFunc(StepPhysics);
    glEnable(GL_DEPTH_TEST);

    glutMainLoop();
    return 0;
}

 实现效果



 

 

  • 大小: 6.4 KB
1
0
分享到:
评论

相关推荐

    质点弹簧模型

    质点弹簧模型是一种常用于模拟物理现象,特别是模拟柔软物体如布料、绳索或弹性材料行为的计算方法。在计算机图形学中,这个模型被广泛应用于游戏开发、动画制作和虚拟现实等领域,以实现逼真的视觉效果。 质点弹簧...

    质点弹簧受重力影响下落的3D效果

    质点弹簧模型是一种常用于计算机图形学中的物理模拟技术,特别是在模拟布料、流体、弹性物体等场景中。在“质点弹簧受重力影响下落的3D效果”这个项目中,我们主要探讨的是如何利用这种算法来创建一个三维环境中,...

    使用弹簧质子模型制作的布料仿真程序.zip

    下面将详细介绍弹簧质点模型及其在布料仿真实现中的应用。 弹簧质点模型是基于物理的模拟系统,其核心思想是将物体视为由许多相互连接的质点组成,这些质点之间通过弹簧进行连接,模拟物体的弹性。每个质点代表物体...

    北航研究生虚拟现实技术课大作业:基于弹簧-质点系统的安全网运动实时3D虚拟仿真

    2. **弹簧-质点系统**:这是一种常见的物理模型,用于描述物体之间的弹性连接。质点代表小且可忽略形状的物体,弹簧则代表连接这些质点的弹性力。在虚拟仿真中,该系统可以模拟物体的动态行为,例如振动、碰撞等。 ...

    opengl绘制弹簧振动(包含全频抗锯齿)

    文件“anti_aliasing”可能包含的是实现这一功能的具体代码示例,可能涵盖了设置OpenGL上下文、创建FBO、配置多重采样参数、绘制弹簧模型以及进行反走样处理的步骤。通过分析和学习这个文件,开发者可以了解到如何将...

    一种基于Vega和OpenGL的的柔性物体建模方法与实现

    通过对质点-弹簧模型的应用以及数学模型的建立,实现了柔性物体的真实动态模拟。这种方法不仅能够提高仿真结果的真实感,还能够有效地应用于虚拟现实、计算机动画等多个领域。通过案例分析表明,该方法能够有效地...

    质点弹簧系统模拟绳子运动.zip

    在质点弹簧模型中,我们可以用弹簧来模拟这种张力。 每个弹簧通常具有两个属性:弹性系数(k)和松弛长度(l0)。弹性系数决定了弹簧恢复原状的力度,而松弛长度是弹簧不受力时的自然长度。当两个相连的质点间距...

    基于改进质点一弹簧模型的图像变形仿真方法 (2009年)

    为了满足仿真系统中对仿真的实时性和视觉效果上的真实性要求,提出了一种改进的质点一弹簧模型,采用...在OpenGL环境下对改进的3D质点一弹簧模型进行了仿真实验,实验结果表明所提模型能逼真地模拟出图像变形的效果。

    基于现代OpenGL的简单物理运动模拟

    根据Nehe第39节教程,使用Opengl核心模式下编程实现匀速,平抛,以及弹簧质点运动。代码根据learnopengl的文件配置方式配置,按着来就好了

    Android 3D游戏开发技术宝典-OpenGL ES 2.0 (吴亚峰) 源代码

    16.4 弹簧质点模型 455 16.4.1 案例运行效果与基本原理 455 16.4.2 具体开发步骤 457 16.5 本章小结 462 第17章 游戏的心脏——物理引擎 463 17.1 物理引擎很重要 463 17.1.1 什么是物理引擎 ...

    北航研究生虚拟现实技术课大作业:基于弹簧-质点系统的窗帘运动实时3D虚拟仿真

    这是我的虚拟现实技术课大作业,用OpenGL和VS2005做的,glut库的各种基本调用都涉及到了,里面有很详细的原理说明。对初学OpenGL者是个很好的例子。 PS: 北航的同学还是乖乖自己写论文吧,毕竟每年真正做出来演示...

    布料运动模拟(质量弹力模型)

    3. **算法实现**:VC++实现的代码需要高效地处理质量弹簧模型的更新、碰撞检测和响应等算法。 五、MosegaardsClothTutorial "**MosegaardsClothTutorial**"可能是一个具体的教程项目,包含源代码、资源文件和说明...

    三维图形程序设计实验指导书.doc

    - **布料模拟**:基于物理的模拟,使用弹簧-质点模型来模拟布料的柔韧性和垂坠感。 - **曲面建模**:包括NURBS(非均匀有理B样条)、Bezier曲面等,用于创建复杂的几何形状。 - **骨架动画**:通过定义骨骼和关节...

    bullet cloth

    为了实现这一目标,开发者们通常会使用一种称为“质点-弹簧系统”的模拟方法,这依赖于物理引擎来计算布料上每个质点间的相互作用力和运动状态。Bullet物理引擎是目前游戏开发中最常用的开源物理SDK之一,由AMD的...

    基于图形硬件的显式织物模拟.pdf

    显式织物模拟通常采用弹簧-质点模型,其中质点代表织物的最小单元,弹簧则模拟纤维间的相互作用力。这种方法简单直观,易于并行化处理。文章中,研究人员利用GPU的并行处理能力,设计了一种新的并行算法。他们将每个...

    皮鞋鞋样CAD系统基础研究与开发.pdf

    这通过建立弹簧质点系统模型,基于能量最小化原理来展开三维曲面到二维平面,从而避免了传统方法中可能出现的鞋面扭曲和重叠问题。 三维编辑技术: 三维编辑技术对于设计人员来说是一个非常重要的工具。它允许设计...

    基于韦尔莱算法的可撕扯幕布模拟

    1. **物理模拟基础**:在模拟幕布的过程中,首先需要建立一个物理模型,包括质点、弹簧和连接它们的约束。质点代表幕布上的每个小块,弹簧则模拟布料间的相互作用力,如拉伸力和剪切力。约束条件确保了布料在空间中...

    实验2_参数关键帧动画设计_visualc++_VS2010_sharpx7u_

    "Mass_Spring-2007.11.23.rar"文件可能是一个关于质点弹簧系统的示例,这是计算机动画中的经典问题。在这个例子中,关键帧动画可以用来模拟物体的弹性运动,通过设置不同时间点的质点位置作为关键帧,然后用插值算法...

    用于软组织变形仿真的GPU辅助能量异步扩散并行计算模型

    3. **异步扩散算法**:通过将质点的机械能引入到异步扩散算法中,能够在每个时间步骤内更精确地计算相邻质点之间的相互作用力,从而实现更加真实的软组织形变效果。 4. **OpenCL加速**:为了满足实时手术仿真对于...

    振动C语言程序代码

    - 简谐振动:以一个质点在弹簧力的作用下进行的往复运动为例,其位移与时间的关系遵循胡克定律和简单的初值条件。 - 阻尼振动:当振动系统存在阻力时,振动会逐渐减小直至停止,其动力学特性需要考虑阻尼系数。 -...

Global site tag (gtag.js) - Google Analytics