[菜鸟每天来段CUDA_C]CUDA实现简单热传导动态模拟

原创 2013年12月03日 15:14:09

本文利用CUDA实现简单二维平面上的热传导模拟。假设有一个矩形房间,将其分成一个网格,在网格中随机散布一些热源,热源有不同的固定温度,然后就可以计算网格中每个单元的温度随时间的变化情况。本文将热传导模型做如下简化:某单元的新时刻温度等于原有温度加其与邻接单元的温差,邻接单元取上下左右四个单元。

温度更新算法:

Step1: 给定一个包含初始输入温度的网格,将其作为热源温度值复制到网格对应的单元;

Step2: 规定一个输入温度网格,根据热传导模型更新温度到输出网格;

Step3: 交换输入温度网格和输出温度网格的数据,转step1.


实验结果如下:下图显示按先后顺序的四个不同时刻的热量分布及变化情况。



主要程序代码如下:

/********************************************************************
*  heatTransfer.cu
*********************************************************************/

#include <stdio.h>
#include <stdlib.h>
#include <cuda_runtime.h>
#include <cutil_inline.h>
#include "CPUAnim.h"

#define DIM 1024
#define MAX_TEMP 1.0f
#define MIN_TEMP 0.0001f
#define SPEED 0.25f

struct DataBlock 
{
	unsigned char	*outputBitmap;
	float			*dev_inSrc;
	float			*dev_outSrc;
	float			*dev_constSrc;
	CPUAnimBitmap	*bitmap;

	cudaEvent_t		start, stop;
	float			totalTime;
	float			frames;
};

template< typename T >
void swap( T& a, T& b ) {
	T t = a;
	a = b;
	b = t;
}


/************************************************************************/
/* Init CUDA                                                            */
/************************************************************************/
bool InitCUDA(void)
{
    ......
}

/************************************************************************/
__global__ void copyConstKernel(float *iPtr, const float *cPtr)
{
	int x = threadIdx.x + blockIdx.x + blockDim.x;
	int y = threadIdx.y + blockIdx.y * blockDim.y;
	int offset = x + y * blockDim.x * gridDim.x;

	if (cPtr[offset] != 0)
	{
		iPtr[offset] = cPtr[offset];
	}
}

__global__ void blendKernel(float* outSrc, const float* inSrc)
{
	int x = threadIdx.x + blockIdx.x * blockDim.x;
	int y = threadIdx.y + blockIdx.y * blockDim.y;
	int offset = x + y * blockDim.x * gridDim.x;

	int left = offset - 1;
	int right = offset + 1;
	if (x == 0) 
	{
		left++; 
	}
	if (x == DIM - 1)
	{
		right--;
	}

	int top = offset - DIM;
	int bottom = offset + DIM;
	if (y == 0)
	{
		top += DIM;
	}
	if (y == DIM -1)
	{
		bottom -= DIM;
	}

	outSrc[offset] = inSrc[offset] + SPEED * (inSrc[top] + inSrc[bottom] 
											+ inSrc[left] + inSrc[right] - 4 * inSrc[offset]);
}

__device__ unsigned char value( float n1, float n2, int hue )
{
      ......
}
__global__ void floatToColor( unsigned char *optr,
							   const float *outSrc ) 
{
    ......
}

void animGPU(DataBlock *d, int ticks)
{
	cutilSafeCall(cudaEventRecord(d->start, 0));
	dim3	blocks(DIM/16, DIM/16);
	dim3	threads(16, 16);
	CPUAnimBitmap	*bitmap = d->bitmap;

	for (int i=0; i<90; i++)
	{
		copyConstKernel<<<blocks, threads>>>(d->dev_inSrc, d->dev_constSrc);
		blendKernel<<<blocks, threads>>>(d->dev_outSrc, d->dev_inSrc);
		swap(d->dev_inSrc, d->dev_outSrc);
	}
	floatToColor<<<blocks, threads>>>(d->outputBitmap, d->dev_inSrc);

	cudaMemcpy(bitmap->get_ptr(), d->outputBitmap, bitmap->image_size(), cudaMemcpyDeviceToHost);
	cutilSafeCall(cudaEventRecord(d->stop, 0));
	cutilSafeCall(cudaEventSynchronize(d->stop));

	float elapsedTime;
	cutilSafeCall(cudaEventElapsedTime(&elapsedTime, d->start, d->stop));
	d->totalTime += elapsedTime;
	++d->frames;

	printf("Average Time per frame: %3.1f ms\n", d->totalTime/d->frames);
}


void animExit(DataBlock *d)
{
	cudaFree(d->dev_inSrc);
	cudaFree(d->dev_outSrc);
	cudaFree(d->dev_constSrc);

	cutilSafeCall(cudaEventDestroy(d->start));
	cutilSafeCall(cudaEventDestroy(d->stop));
}
/************************************************************************/

int main(int argc, char* argv[])
{

	if(!InitCUDA()) {
		return 0;
	}

	DataBlock	data;
	CPUAnimBitmap bitmap(DIM, DIM, &data);
	data.bitmap = &bitmap;
	data.totalTime = 0;
	data.frames = 0;

	cutilSafeCall(cudaEventCreate(&data.start));
	cutilSafeCall(cudaEventCreate(&data.stop));

	cutilSafeCall(cudaMalloc((void**)&data.outputBitmap, data.bitmap->image_size()));

	cutilSafeCall(cudaMalloc((void**)&data.dev_inSrc, data.bitmap->image_size()));
	cutilSafeCall(cudaMalloc((void**)&data.dev_outSrc, data.bitmap->image_size()));
	cutilSafeCall(cudaMalloc((void**)&data.dev_constSrc, data.bitmap->image_size()));

	float *temp = (float*)malloc(bitmap.image_size());
	for (int i=0; i<DIM*DIM; i++)
	{
		temp[i] = 0;
		int x = i % DIM;
		int y = i / DIM;

		if ((x>300) && (x<600) && (y>301) && (y<601))
		{
			temp[i] = MAX_TEMP;
		}
	}
	temp[DIM*100+100] = (MAX_TEMP + MIN_TEMP) / 2;
    temp[DIM*700+100] = MIN_TEMP;
	temp[DIM*300+300] = MIN_TEMP;
	temp[DIM*200+700] = MIN_TEMP;

	for (int y=800; y<900; y++)
	{
		for (int x=400; x<500; x++)
		{
			temp[x+y*DIM] = MIN_TEMP;
		}
	}

	cutilSafeCall(cudaMemcpy(data.dev_constSrc, temp, bitmap.image_size(), cudaMemcpyHostToDevice));
	
	for(int y=800; y<DIM; y++)
	{
		for (int x=0; x<200; x++)
		{
			temp[x+y*DIM] = MAX_TEMP;
		}
	}
	cutilSafeCall(cudaMemcpy(data.dev_inSrc, temp, bitmap.image_size(), cudaMemcpyHostToDevice));
	free(temp);
	bitmap.anim_and_exit((void(*)(void*, int))animGPU, (void(*)(void*))animExit);
	return 0;
}


[菜鸟每天来段CUDA_C]CUDA实现向量的点积运算

本文利用CUDA实现向量的点积运算。主要思想是:每个线程计算两个相应元素的乘积,然后利用共享内存(__shared__)累计每个线程的结果,得到一个线程块内向量的部分内积,最后在CPU上对个线程块的结...
  • jonny_super
  • jonny_super
  • 2013年11月22日 14:07
  • 1789

FDM之二维静态热传导--含有不同传导系数K

% Geol 5200, Humphrey 2013 % 2D, steady state temperature.  Simplest approach % full (not sparse) ma...
  • conjimmy
  • conjimmy
  • 2015年03月11日 07:15
  • 445

小白福利!思维导图超详细绘制教程

思维导图,也就是心智图或者脑图,是近年来十分流行的一种思维管理方法。有的人喜欢把图画的花花绿绿的,便于自己加深记忆,这并没有错,但是思维导图本来就是用来提高效率的工具,一张这么“漂亮”的图画下来,时间...
  • mastermindmap
  • mastermindmap
  • 2017年11月27日 10:41
  • 109

二维热传导方程求解

  • 2017年09月21日 23:24
  • 3.16MB
  • 下载

[菜鸟每天来段CUDA_C]CppIntegration在C++程序中引用CUDA程序

本文主要实现在C++程序中引用CUDA程序,主要意义是使顺序定义的数据能在CUDA程序中并行执行,然后返回结果。 程序主要包括main.cpp                           ...
  • jonny_super
  • jonny_super
  • 2013年11月19日 21:38
  • 1144

[菜鸟每天来段CUDA_C]CUDA与OpenGL互操作

本文要解决的问题是如何实现CUDA和OpenGL的互操作,使得GPU能够将通用计算的运算结果交给OpenGL进行绘制。 本文的应用程序主要包括两个方面: 1.      使用CUDA核函数生成图像...
  • jonny_super
  • jonny_super
  • 2013年12月08日 15:28
  • 2254

思维导图

思维导图是一个描述和引导思维的工具。有关思维导图请参看维基百科“思维导图”,“心智图”词条。 绘制思维导图的工具非常多(百科中有介绍),笔者使用Xmind。 网盘下...
  • yanwushu
  • yanwushu
  • 2012年07月30日 11:38
  • 1687

Java技术导图

简单的画了一张有关Java的技术导图,其中不免纰漏,作为自己学习的一个指导方向。...
  • luohuacanyue
  • luohuacanyue
  • 2012年12月08日 15:12
  • 1596

基于MATLAB的有限差分法解决二位瞬态导热问题

本程序解决的是应对平板的二维导热问题,可以定义边界温度及平板材料,获得一定时间后的温度云图。它只考虑热传导,没有考虑对流传热及热辐射,自己做个记录吧。 程序代码如下: clc; clea...
  • sealingj
  • sealingj
  • 2017年01月08日 10:05
  • 539

实现波纹效果

cuda用来算像素颜色+openGL用来显示(mdzz) 波纹的颜色函数直接用的书上建议不要花太多时间去理解的函数 放代码#include "cuda_runtime.h" #include "d...
  • u011229924
  • u011229924
  • 2017年04月12日 18:24
  • 440
内容举报
返回顶部
收藏助手
不良信息举报
您举报文章:[菜鸟每天来段CUDA_C]CUDA实现简单热传导动态模拟
举报原因:
原因补充:

(最多只允许输入30个字)