梯形积分法【OpenMP实现】多个版本

简述

会用多个版本来写。操作是用命令行来操作。
除了第二版本是在内容上基于第一个版本的完善，其他都是降低算法的复杂程度的，可以放心阅读

文章目录

简述
版本1
版本2
版本3
版本4

版本1

要求n被thread_count整除 不然可以会有计算错误。

#include <iostream>
#include <omp.h>#define FUN(x) (x * x)using namespace std;
#pragma warning(disable : 4996)
void Trap(double a, double b, int n, double *global_result_p);
int main(int argc, char **argv) {if (argc < 5) return 0;double global_result = 0.0;double a, b;int n;int thread_count = strtol(argv[1], NULL, 10);a = strtod(argv[2], NULL);b = strtod(argv[3], NULL);n = strtol(argv[4], NULL, 10);// n 必须被 thread_count 整除
#pragma omp parallel num_threads (thread_count)Trap(a, b, n, &global_result);cout << "With n = " << n << " trapzoids, our estimate\n";cout << "of the integral from " << a << " to " << b <<" = "<< global_result<< endl;
}void Trap(double a, double b, int n, double *global_result_p) {double h, x, my_result;double local_a, local_b;int i, local_n;int my_rank = omp_get_thread_num();int thread_count = omp_get_num_threads();h = (b - a) / n;local_n = n / thread_count;  // same in all threadslocal_a = a + my_rank * local_n * h;local_b = local_a + local_n * h;x = local_a; // initial xmy_result = (FUN(local_a) + FUN(local_b)) / 2;for (i = 1; i < local_n; ++i) {x += h;my_result += FUN(x);}my_result *= h;
# pragma omp critical *global_result_p += my_result;
}

编译好源文件之后，就用下面的命令执行

这里是用4个线程，来计算在0，1区间上的x2x^2x2的用n=64划分的梯形积分公式。结果接近1/3 可以通过提高n来提高精度。但这个版本只能计算n为thread_count的整数倍的时候

PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 64
With n = 64 trapzoids, our estimate
of the integral from 0 to 1 = 0.333374

版本2

可以使用不被thread_count整除的n
改进代码： 添加了一个变量left表示剩余。剩余的让前left个线程每个再多算点。（一般来说线程数目都会比n要小很多。（不然在算法上就会不精确（n不够大的话）。））

left = n % thread_count;
if (my_rank < left) local_n += 1;
if (my_rank < left){ local_a = a + my_rank * local_n * h; }
else { local_a = a + left * (local_n + 1) * h + (my_rank - left) * local_n * h; }

实际的代码:

#include <iostream>
#include <omp.h>#define FUN(x) (x * x)using namespace std;
#pragma warning(disable : 4996)
void Trap(double a, double b, int n, double *global_result_p);
int main(int argc, char **argv) {if (argc < 5) return 0;double global_result = 0.0;double a, b;int n;int thread_count = strtol(argv[1], NULL, 10);a = strtod(argv[2], NULL);b = strtod(argv[3], NULL);n = strtol(argv[4], NULL, 10);#pragma omp parallel num_threads (thread_count)Trap(a, b, n, &global_result);cout << "With n = " << n << " trapzoids, our estimate\n";cout << "of the integral from " << a << " to " << b <<" = "<< global_result<< endl;
}void Trap(double a, double b, int n, double *global_result_p) {double h, x, my_result;double local_a, local_b;int i, local_n, left;int my_rank = omp_get_thread_num();int thread_count = omp_get_num_threads();h = (b - a) / n;local_n = n / thread_count;  // same in all threadsleft = n % thread_count;if (my_rank < left) local_n += 1;if (my_rank < left){ local_a = a + my_rank * local_n * h; }else { local_a = a + left * (local_n + 1) * h + (my_rank - left) * local_n * h; }local_b = local_a + local_n * h;x = local_a; // initial xmy_result = (FUN(local_a) + FUN(local_b)) / 2;for (i = 1; i < local_n; ++i) {x += h;my_result += FUN(x);}my_result *= h;
# pragma omp critical *global_result_p += my_result;
}

发现精度不断的提升，在只改变n的数目的情况下。

PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 64
With n = 64 trapzoids, our estimate
of the integral from 0 to 1 = 0.333374
PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 65
With n = 65 trapzoids, our estimate
of the integral from 0 to 1 = 0.333373
PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 66
With n = 66 trapzoids, our estimate
of the integral from 0 to 1 = 0.333372
PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 67
With n = 67 trapzoids, our estimate
of the integral from 0 to 1 = 0.33337

版本3

对于部分对指针不太熟的人来说，我们可以做出下面的改进
修改部分：

#pragma omp parallel num_threads (thread_count){double my_result = Trap(a, b, n);
#pragma omp critical global_result += my_result;}

#include <iostream>
#include <omp.h>#define FUN(x) (x * x)using namespace std;
#pragma warning(disable : 4996)
double Trap(double a, double b, int n);
int main(int argc, char **argv) {if (argc < 5) return 0;double global_result = 0.0;double a, b;int n;int thread_count = strtol(argv[1], NULL, 10);a = strtod(argv[2], NULL);b = strtod(argv[3], NULL);n = strtol(argv[4], NULL, 10);#pragma omp parallel num_threads (thread_count){double my_result = Trap(a, b, n);
#pragma omp critical global_result += my_result;}cout << "With n = " << n << " trapzoids, our estimate\n";cout << "of the integral from " << a << " to " << b <<" = "<< global_result<< endl;
}double Trap(double a, double b, int n) {double h, x, my_result;double local_a, local_b;int i, local_n, left;int my_rank = omp_get_thread_num();int thread_count = omp_get_num_threads();h = (b - a) / n;local_n = n / thread_count;  // same in all threadsleft = n % thread_count;if (my_rank < left) local_n += 1;if (my_rank < left){ local_a = a + my_rank * local_n * h; }else { local_a = a + left * (local_n + 1) * h + (my_rank - left) * local_n * h; }local_b = local_a + local_n * h;x = local_a; // initial xmy_result = (FUN(local_a) + FUN(local_b)) / 2;for (i = 1; i < local_n; ++i) {x += h;my_result += FUN(x);}my_result *= h;return my_result;
}

结果

PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 67
With n = 67 trapzoids, our estimate
of the integral from 0 to 1 = 0.33337
PS D:\C++\VS\repo\OpenMP-TEST\Debug> ./OpenMP-TEST 4 0 1 68
With n = 68 trapzoids, our estimate
of the integral from 0 to 1 = 0.333369
PS D:\C++\VS\repo\OpenMP-TEST\Debug>

版本4

使用规约变量和规约操作
会让代码更加简单

#pragma omp parallel num_threads (thread_count) reduction(+:global_result){global_result += Trap(a, b, n);}

跟版本3是一样的。就是让代码更加简单了而已，用了map-reduce的思路
reduction(<operator>: <variable list>) 前面的放操作符，后面的设置变量。
但是有需要注意的，如果是减法就需要小心了。在每个线程上确实是使用了减法，但是全局上也是可能做了减法，(例如 -1 - (-5)= 4）这样的事情是有可能发生的。一般来说，满足交换律的都不用担心这点hh

梯形积分法【OpenMP实现】多个版本相关推荐

关于mpi的理论知识以及编写程序来实现数据积分中的梯形积分法。
几乎所有人的第一个程序是从"hello,world"程序开始学习的 #include "mpi.h" #include <stdio.h> int ...
MPI——实现梯形积分法
用MPI来实现梯度积分法梯度积分法我们可以用积分法来估计函数y=f(x)y=f(x)y=f(x)的图像中,两条垂直线与xxx轴之间的区间大小基本思想:将xxx轴上的区间划分为n个等长的子区间.然 ...
复化梯形公式matlab程序_【HPC】高性能数值计算-梯形面积法
原理将 x 轴上的区间划分为 n 个等长的子区间,有 , , ,......, , 则梯形面积可表示为首先很容易可以得出串行计算的代码: /*Input a, b, n*/ h = (b - a) ...
并行程序设计方法实验（包括openmp、向量化实现pi计算、SPECOMP2012测试、矩阵乘法优化）
目录一.实验环境二.专题一之积分计算圆周率 2.1向量优化 2.2 OpenMP优化三.专题二之测试SPECOMP2012 3.1初步了解SPECOMP 3.2系统基本配置 3.3实践 3.3. ...
Ubuntu19.10下TAU的配置及梯形积分法的实现
目录前言一. 下载并安装TAU的前置事项二. 编译并安装TAU 三. 编写代码并运行 (1)串行部分 (2)并行部分 (3)结果分析前言 TAU是一种可以在Ubuntu下对并行运算的进程进行性 ...
c语言用蒙特卡洛算法求定积分,C++计算定积分、不定积分、蒙特卡洛积分法
#include"integral.h" double integral::definiteIntegral(double(*f)(double), double a, doubl ...
并行程序设计——OMP编程
并行程序设计--OMP编程实验一实验内容分别实现课件中的梯形积分法的Pthread.OpenMP版本,熟悉并掌握OpenMP编程方法,探讨两种编程方式的异同. 实验代码 OpenMP编程 #in ...
mpi并行 java_【并行计算】用MPI进行分布式内存编程（一）
通过上一篇关于并行计算准备部分的介绍,我们知道MPI(Message-Passing-Interface 消息传递接口)实现并行是进程级别的,通过通信在进程之间进行消息传递.MPI并不是一种新的开发语 ...
matlab符号值,matlab 符号积分和数值积分
f=(x)sin(x)+cos(x);定义被积函数Intf=quad(f,0,1)%对f进行积分,下限为0,上限为1www.mh456.com防采集. 其实,你仔细观察一下会发现,用第二种方法得到百x ...

梯形积分法【OpenMP实现】多个版本

简述

文章目录

版本1

版本2

版本3

版本4

梯形积分法【OpenMP实现】多个版本相关推荐

最新文章

热门文章