文章目录
- 1. 学习目的
- 2. 阶段成果
- 2.1 NVVP 性能探查
- 2.2 测试编译环境
- 2.3 测试样例
- 3 学习过程
- 3.1 提问DeepSeek
- 3.2 最终代码
- 4. 体会
1. 学习目的
最近在学习cuda，准备给我的taskBus SDR添加CUDA的模块支持，以便可以用PC机压榨山寨 B210那56M的带宽。
因为年纪比较大了，所以不打算看书了，直接DeepSeek安排上，结果发现学习效率杠杠的！
2. 阶段成果
我用自己的文物显卡 NVIDIA GeForce GTX 750 Ti (Compute Capability 5.0)，只有2G的显存，1G的主频，做时域512点fir卷积，竟然只需0.65秒就跑完100M float数据（411MB文件大小）。
关键技术是使用了三个线程进行读写处理分离。
2.1 NVVP 性能探查
使用NVVP，可以看到GPU的占空比还是很紧凑。
主要时间花在了fir_kernel，这是我们希望的。同时，部分时间在进行GPU内部的memCpy，这也是值得进一步优化的方向。
由于磁盘IO完全和CPU计算独立，所以两波计算之间几乎没有空闲。
2.2 测试编译环境
Manjaro 20250321，CUDA12
Intel Core i7 6500K @ 4.0GHz
编译命令：
nvcc -arch=sm_50 cuda_filter.cu -o cuda_filter
注意要使用 -arch=sm_50 选择老的计算兼容选项，否则CUDA启动后会报错。
2.3 测试样例
测试样例是用 AudaCity录制的音频，411MB大小，100MB左右的float32样点（单声道）
滤波器对数据（轨道1）做512点的平均（192KHz采样率下，大概是2秒多）。
3 学习过程
3.1 提问DeepSeek
提问：
使用C语言调用CUDA实现一个连续FIR滤波算法。用户使用该算法的步骤是：
- 打开文件 test.pcm。
- 循环每次读取一定个数（N个）float值，直到文件尾部。
- 对每次读取的float值，调用 push_and_filter() 方法。该方法是你需要着重考虑设计、实现的，作用是接着上次的状态，继续执行连续滤波。滤波器是一个提前定义好的数组，有M个元素，float类型，M目前设置合理范围是16-2048。
- 每次滤波的输出，存储在 out.pcm,连续写入。
- 滤波缓存一开始为全0,尾部也全部清零后结束。
请根据CUDA的特性，分析并设计每次读取的N的个数，让计算尽可能快，并给出全部代码。
回答：
#include <stdio.h>
#include <stdlib.h>
#include <cuda_runtime.h>#define MAX_FILTER_LENGTH 2048
#define OPTIMAL_N 4096 // 根据GPU架构优化的读取块大小
#define BLOCK_SIZE 256 // CUDA线程块大小typedef struct {float* buffer;int size;
} CircularBuffer;__constant__ float d_filter[MAX_FILTER_LENGTH]; // 滤波器系数存储在常量内存// CUDA核函数，使用共享内存优化
__global__ void fir_kernel(const float* input, float* output, int N, int M) {extern __shared__ float s_window[];const int tid = threadIdx.x;const int bid = blockIdx.x;const int gid = bid * blockDim.x + tid;// 每个块加载自己需要的输入窗口const int window_size = blockDim.x + MAX_FILTER_LENGTH - 1;const int input_offset = bid * blockDim.x;// 协作加载数据到共享内存for(int i = tid; i < window_size; i += blockDim.x) {const int load_pos = input_offset + i;s_window[i] = (load_pos < N + M - 1) ? input[load_pos] : 0.0f;}__syncthreads();if(gid < N) { // 确保不越界float sum = 0.0f;for(int j = 0; j < M; j++) {sum += s_window[tid + j] * d_filter[j];}output[gid] = sum;}
}void push_and_filter(CircularBuffer* state, const float* input, int input_len,float** output, int* output_len, int M, cudaStream_t stream) {// 准备扩展输入缓冲区（历史数据 + 新输入）const int total_len = state->size + input_len;float* extended_input = (float*)malloc(total_len * sizeof(float));// 合并历史数据和当前输入memcpy(extended_input, state->buffer, state->size * sizeof(float));memcpy(extended_input + state->size, input, input_len * sizeof(float));// 准备设备内存float *d_input, *d_output;const int output_size = input_len;cudaMalloc(&d_input, total_len * sizeof(float));cudaMalloc(&d_output, output_size * sizeof(float));// 异步内存拷贝和计算cudaMemcpyAsync(d_input, extended_input, total_len * sizeof(float),cudaMemcpyHostToDevice, stream);const dim3 grid_dim((output_size + BLOCK_SIZE - 1) / BLOCK_SIZE);const dim3 block_dim(BLOCK_SIZE);const size_t shared_mem = (BLOCK_SIZE + MAX_FILTER_LENGTH - 1) * sizeof(float);fir_kernel<<<grid_dim, block_dim, shared_mem, stream>>>(d_input, d_output, output_size, M);// 回收历史数据（保留最后M-1个样本）free(state->buffer);state->size = M - 1;state->buffer = (float*)malloc(state->size * sizeof(float));memcpy(state->buffer, extended_input + input_len, state->size * sizeof(float));// 异步获取结果*output = (float*)malloc(output_size * sizeof(float));cudaMemcpyAsync(*output, d_output, output_size * sizeof(float),cudaMemcpyDeviceToHost, stream);cudaStreamSynchronize(stream);*output_len = output_size;free(extended_input);cudaFree(d_input);cudaFree(d_output);
}int main() {// 初始化滤波器系数（示例值）const int M = 512;float h_filter[M];for(int i = 0; i < M; i++) h_filter[i] = 1.0f/M;cudaMemcpyToSymbol(d_filter, h_filter, M * sizeof(float));// 创建CUDA流和状态缓冲区cudaStream_t stream;cudaStreamCreate(&stream);CircularBuffer state = {NULL, 0};FILE *fin = fopen("test.pcm", "rb");FILE *fout = fopen("out.pcm", "wb");float *input_buf = (float*)malloc(OPTIMAL_N * sizeof(float));while(1) {const size_t read_cnt = fread(input_buf, sizeof(float), OPTIMAL_N, fin);if(read_cnt == 0) break;float* output_buf;int output_len;push_and_filter(&state, input_buf, read_cnt, &output_buf, &output_len, M, stream);fwrite(output_buf, sizeof(float), output_len, fout);free(output_buf);} // 清理资源free(input_buf);free(flush_buf);free(state.buffer);fclose(fin);fclose(fout);cudaStreamDestroy(stream);return 0;
}
而后，利用多线程修改，多次互动后，得到最终代码：
3.2 最终代码
#include <stdio.h>
#include <stdlib.h>
#include <cuda_runtime.h>
#include <atomic>
#include <chrono>
#include <thread>
#define MAX_FILTER_LENGTH 2048
#define OPTIMAL_N 1024*1024 // Optimal read block size based on GPU architecture
#define BLOCK_SIZE 256 // CUDA thread block size
#define NUM_STREAMS 4
#define CHECK_CUDA(call) { \cudaError_t err = call; \if (err != cudaSuccess) { \fprintf(stderr, "CUDA error at %s:%d: %s\n", __FILE__, __LINE__, cudaGetErrorString(err)); \exit(EXIT_FAILURE); \} \
}typedef struct {float* buffer;int size;int max_size;
} CircularBuffer;__constant__ float d_filter[MAX_FILTER_LENGTH]; // Filter coefficients in constant memory// CUDA kernel with shared memory optimization
__global__ void fir_kernel(const float* input, float* output, int N, int M) {extern __shared__ float s_window[];const int tid = threadIdx.x;const int bid = blockIdx.x;const int gid = bid * blockDim.x + tid;// Each block loads its needed input windowconst int window_size = blockDim.x + MAX_FILTER_LENGTH - 1;const int input_offset = bid * blockDim.x;// Cooperative data loading into shared memoryfor(int i = tid; i < window_size; i += blockDim.x) {const int load_pos = input_offset + i;s_window[i] = (load_pos < N + M - 1) ? input[load_pos] : 0.0f;}__syncthreads();if(gid < N) { // Ensure we don't go out of boundsfloat sum = 0.0f;for(int j = 0; j < M; j++) {sum += s_window[tid + j] * d_filter[j];}output[gid] = sum;}
}
//用于拼接前序状态的输入缓存
int mem_ext_input_size = 0;
float* mem_ext_input = 0;
//用于CUDA设备侧的输入输出缓存
float *d_input = 0, *d_output = 0;
int d_inputlen = 0, d_outputlen = 0;std::atomic<bool> read_finished(false),deal_finished(false);
std::atomic<long long> read_pos(0);
std::atomic<long long> deal_pos(0);
std::atomic<long long> write_pos(0);//连续滤波函数
void push_and_filter(CircularBuffer* state, const float* input, int input_len,float* output, int * outlen, int M, cudaStream_t stream) {// Prepare extended input buffer (history data + new input)int ext_input_len = state->size + input_len;if (mem_ext_input_size < ext_input_len) {if (mem_ext_input)free(mem_ext_input);mem_ext_input = (float*)malloc(ext_input_len * sizeof(float));mem_ext_input_size = ext_input_len;}// Combine history data with current inputmemcpy(mem_ext_input, state->buffer, state->size * sizeof(float));memcpy(mem_ext_input + state->size, input, input_len * sizeof(float));// Prepare device memoryconst int output_size = input_len;if (d_inputlen < ext_input_len) {if (d_input)cudaFree(d_input);CHECK_CUDA(cudaMalloc(&d_input, ext_input_len * sizeof(float)));d_inputlen = ext_input_len;}if (d_outputlen < output_size) {if (d_output)cudaFree(d_output);CHECK_CUDA(cudaMalloc(&d_output, output_size * sizeof(float)));d_outputlen = output_size;}// Synchronous memory copy//CHECK_CUDA(cudaMemcpy(d_input, mem_ext_input, ext_input_len * sizeof(float), cudaMemcpyHostToDevice));// Async memory operations with streamCHECK_CUDA(cudaMemcpyAsync(d_input, mem_ext_input, ext_input_len * sizeof(float), cudaMemcpyHostToDevice, stream));const dim3 grid_dim((output_size + BLOCK_SIZE - 1) / BLOCK_SIZE);const dim3 block_dim(BLOCK_SIZE);const size_t shared_mem = (BLOCK_SIZE + MAX_FILTER_LENGTH - 1) * sizeof(float);// Kernel launch with default streamfir_kernel<<<grid_dim, block_dim, shared_mem ,stream>>>(d_input, d_output, output_size, M);// Recycle history data (keep last M-1 samples)state->size = M - 1;if (state->max_size < state->size){if (state->buffer )free (state->buffer );state->buffer = (float*)malloc(state->size * sizeof(float));state->max_size = state->size;}memcpy(state->buffer, mem_ext_input + input_len, state->size * sizeof(float));//CHECK_CUDA(cudaMemcpy(output, d_output, output_size * sizeof(float), cudaMemcpyDeviceToHost));CHECK_CUDA(cudaMemcpyAsync(output, d_output, output_size * sizeof(float), cudaMemcpyDeviceToHost, stream));*outlen = output_size;CHECK_CUDA(cudaDeviceSynchronize());return;
}int main() {// Initialize filter coefficients (example values)const int M = 512;float h_filter[M];for(int i = 0; i < M; i++) h_filter[i] = 1.0f/M;CHECK_CUDA(cudaMemcpyToSymbol(d_filter, h_filter, M * sizeof(float)));// Initialize state bufferCircularBuffer state = {NULL, 0 , 0};const int buf_stage = 8;float (*input_buf)[OPTIMAL_N] = new float [buf_stage][OPTIMAL_N];float (*output_buf)[OPTIMAL_N] = new float [buf_stage][OPTIMAL_N];size_t (*input_len) = new size_t [buf_stage];int (*output_len) = new int [buf_stage];//Read threadstd::thread read_thread([&]()->void{FILE *fin = fopen("test.pcm", "rb");while(!read_finished) {if (deal_pos + buf_stage - 2 <= read_pos){//Sleepstd::chrono::milliseconds dura(1);std::this_thread::sleep_for(dura);continue;}const long long pos = read_pos % buf_stage;//Readinput_len[pos] = fread(input_buf[pos], sizeof(float), OPTIMAL_N, fin);if(input_len[pos] < OPTIMAL_N) read_finished = true;if(++read_pos % 10 ==0) fprintf(stderr,"ReadPos = %lld\n",(long long)read_pos);} fclose(fin);});//write threadstd::thread write_thread([&]()->void{FILE *fout = fopen("out.pcm", "wb");while((!deal_finished) || write_pos < deal_pos) {if (write_pos >= deal_pos){//Sleepstd::chrono::milliseconds dura(1);std::this_thread::sleep_for(dura);continue;}const long long pos = write_pos % buf_stage;if (output_len[pos])fwrite(output_buf[pos], sizeof(float), output_len[pos], fout);if (++write_pos % 10 ==0)fprintf(stderr,"WritePos = %lld\n",(long long)write_pos);} fclose(fout);});cudaStream_t streams[buf_stage];for (int i = 0; i < buf_stage; i++) {cudaStreamCreate(&streams[i]);}//Deal in main Threadwhile ((!read_finished) || deal_pos < read_pos ){if (deal_pos >= read_pos){//Sleepstd::chrono::milliseconds dura(1);std::this_thread::sleep_for(dura);continue;}const long long pos = deal_pos % buf_stage;cudaStream_t stream = streams[pos % NUM_STREAMS];push_and_filter(&state, input_buf[pos], input_len[pos], output_buf[pos], &(output_len[pos]), M,stream); if (++deal_pos % 10==0)fprintf(stderr,"DealPos = %lld\n",(long long)deal_pos);}deal_finished = true;read_thread.join();write_thread.join();fprintf(stderr,"Deal Finished.\n");if (state.buffer) free(state.buffer);if (mem_ext_input) free(mem_ext_input);if (d_input) CHECK_CUDA(cudaFree(d_input));if (d_output) CHECK_CUDA(cudaFree(d_output));delete []input_buf;delete []output_buf;delete []input_len;delete []output_len;for (int i = 0; i < NUM_STREAMS; i++) {CHECK_CUDA(cudaStreamDestroy(streams[i]));}return 0;
}4. 体会
从不熟悉CUDA，到使用CUDA完成了性能还乏善可陈的滤波器，用了4个小时。同时，通过询问DeepSeek如何理解stream流式处理、如何使用Nvidia Visual Profiler的优化指南，了解到很多CUDA相关的周边知识。可见，AI已经使得我这个50岁的人都能学到新东西，更别提年轻人了。现在是只有想不到，没有做不到。
】安全审计与日志管理：为电商平台筑牢安全防线)
